1 #ifndef AMREX_WRITE_BINARY_PARTICLE_DATA_HDF5_H
2 #define AMREX_WRITE_BINARY_PARTICLE_DATA_HDF5_H
7 hid_t attr, attr_space;
10 attr_space = H5Screate_simple(1, &dims, NULL);
12 attr = H5Acreate(loc, name, dtype, attr_space, H5P_DEFAULT, H5P_DEFAULT);
14 amrex::Print() <<
" Error with H5Acreate [" << name <<
"]\n";
18 ret = H5Awrite(attr, dtype, (
void*)data);
20 amrex::Print() <<
" Error with H5Awrite [" << name <<
"]\n";
30 hid_t attr, atype, space;
36 space = H5Screate(H5S_SCALAR);
37 atype = H5Tcopy(H5T_C_S1);
38 H5Tset_size(atype, strlen(str)+1);
39 H5Tset_strpad(atype,H5T_STR_NULLTERM);
40 attr = H5Acreate(loc, name, atype, space, H5P_DEFAULT, H5P_DEFAULT);
42 amrex::Print() <<
" Error with H5Acreate [" << name <<
"]\n";
46 ret = H5Awrite(attr, atype, str);
48 amrex::Print() <<
" Error with H5Awrite [" << name <<
"]\n";
59 static int ReadHDF5Attr(hid_t loc,
const char *name,
void *data, hid_t dtype)
64 attr = H5Aopen(loc, name, H5P_DEFAULT);
66 amrex::Print() <<
" Error with H5Aopen [" << name <<
"]\n";
70 ret = H5Aread(attr, dtype, data);
72 amrex::Print() <<
" Error with H5Aread [" << name <<
"]\n";
86 H5Pset_fapl_mpio(fapl, comm, MPI_INFO_NULL);
89 int alignment = 16 * 1024 * 1024;
90 int blocksize = 4 * 1024 * 1024;
91 H5Pset_alignment(fapl, alignment, alignment);
92 H5Pset_meta_block_size(fapl, blocksize);
95 H5Pset_coll_metadata_write(fapl,
true);
96 H5Pset_all_coll_metadata_ops(fapl,
true);
99 H5AC_cache_config_t cache_config;
100 cache_config.version = H5AC__CURR_CACHE_CONFIG_VERSION;
101 H5Pget_mdc_config(fapl, &cache_config);
102 cache_config.set_initial_size = 1;
103 cache_config.initial_size = 16 * 1024 * 1024;
104 cache_config.evictions_enabled = 0;
105 cache_config.incr_mode = H5C_incr__off;
106 cache_config.flash_incr_mode = H5C_flash_incr__off;
107 cache_config.decr_mode = H5C_decr__off;
108 H5Pset_mdc_config (fapl, &cache_config);
110 H5Pset_fapl_sec2(fapl);
114 template <class PC, class F, std::enable_if_t<IsParticleContainer<PC>::value,
int> foo = 0>
116 const std::string& dir,
const std::string& name,
117 const Vector<int>& write_real_comp,
118 const Vector<int>& write_int_comp,
119 const Vector<std::string>& real_comp_names,
120 const Vector<std::string>& int_comp_names,
121 const std::string& compression,
122 F&&
f,
bool is_checkpoint)
127 AMREX_ASSERT(
sizeof(
typename PC::ParticleType::RealType) == 4 ||
128 sizeof(
typename PC::ParticleType::RealType) == 8);
130 constexpr
int NStructReal = PC::NStructReal;
131 constexpr
int NStructInt = PC::NStructInt;
139 #ifdef AMREX_USE_HDF5_ASYNC
142 async_vol_es_wait_particle();
145 es_par_g = H5EScreate();
149 std::string pdir = dir;
150 if ( ! pdir.empty() && pdir[pdir.size()-1] !=
'/') { pdir +=
'/'; }
153 if ( ! pc.GetLevelDirectoriesCreated()) {
165 hid_t fid, grp, fapl, comp_dtype;
171 comp_dtype = H5Tcreate (H5T_COMPOUND, 2 * AMREX_SPACEDIM *
sizeof(
int));
172 if (1 == AMREX_SPACEDIM) {
173 H5Tinsert (comp_dtype,
"lo_i", 0 *
sizeof(
int), H5T_NATIVE_INT);
174 H5Tinsert (comp_dtype,
"hi_i", 1 *
sizeof(
int), H5T_NATIVE_INT);
176 else if (2 == AMREX_SPACEDIM) {
177 H5Tinsert (comp_dtype,
"lo_i", 0 *
sizeof(
int), H5T_NATIVE_INT);
178 H5Tinsert (comp_dtype,
"lo_j", 1 *
sizeof(
int), H5T_NATIVE_INT);
179 H5Tinsert (comp_dtype,
"hi_i", 2 *
sizeof(
int), H5T_NATIVE_INT);
180 H5Tinsert (comp_dtype,
"hi_j", 3 *
sizeof(
int), H5T_NATIVE_INT);
182 else if (3 == AMREX_SPACEDIM) {
183 H5Tinsert (comp_dtype,
"lo_i", 0 *
sizeof(
int), H5T_NATIVE_INT);
184 H5Tinsert (comp_dtype,
"lo_j", 1 *
sizeof(
int), H5T_NATIVE_INT);
185 H5Tinsert (comp_dtype,
"lo_k", 2 *
sizeof(
int), H5T_NATIVE_INT);
186 H5Tinsert (comp_dtype,
"hi_i", 3 *
sizeof(
int), H5T_NATIVE_INT);
187 H5Tinsert (comp_dtype,
"hi_j", 4 *
sizeof(
int), H5T_NATIVE_INT);
188 H5Tinsert (comp_dtype,
"hi_k", 5 *
sizeof(
int), H5T_NATIVE_INT);
192 Vector<std::map<std::pair<int, int>,
typename PC::IntVector > >
193 particle_io_flags(pc.GetParticles().size());
194 for (
int lev = 0; lev < pc.GetParticles().
size(); lev++)
196 const auto& pmap = pc.GetParticles(lev);
197 for (
const auto& kv : pmap)
199 auto& flags = particle_io_flags[lev][kv.first];
206 if(pc.GetUsePrePost())
208 nparticles = pc.GetNParticlesPrePost();
209 maxnextid = pc.GetMaxNextIDPrePost();
214 maxnextid = PC::ParticleType::NextID();
216 PC::ParticleType::NextID(maxnextid);
220 std::string HDF5FileName = pdir;
221 if ( ! HDF5FileName.empty() && HDF5FileName[HDF5FileName.size()-1] !=
'/')
224 HDF5FileName += name;
225 HDF5FileName +=
".h5";
226 pc.HdrFileNamePrePost = HDF5FileName;
233 char setstripe[1024];
234 int stripe_count = 128;
236 char *stripe_count_str = getenv(
"HDF5_STRIPE_COUNT");
237 char *stripe_size_str = getenv(
"HDF5_STRIPE_SIZE");
238 if (stripe_count_str) {
239 stripe_count = atoi(stripe_count_str);
242 if (stripe_size_str) {
243 stripe_size = atoi(stripe_size_str);
246 if (set_stripe == 1) {
247 snprintf(setstripe,
sizeof setstripe,
"lfs setstripe -c %d -S %dm %s", stripe_count, stripe_size, pdir.c_str());
248 std::cout <<
"Setting stripe parameters for HDF5 output: " << setstripe << std::endl;
252 fid = H5Fcreate(HDF5FileName.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
261 std::string versionName = is_checkpoint ? PC::CheckpointVersion() : PC::PlotfileVersion();
262 if (
sizeof(
typename PC::ParticleType::RealType) == 4) {
263 versionName +=
"_single";
265 versionName +=
"_double";
270 int num_output_real = 0;
271 for (
int i = 0; i < pc.NumRealComps() + NStructReal; ++i) {
272 if (write_real_comp[i]) { ++num_output_real; }
275 int num_output_int = 0;
276 for (
int i = 0; i < pc.NumIntComps() + NStructInt; ++i) {
277 if (write_int_comp[i]) { ++num_output_int; }
281 int ndim = AMREX_SPACEDIM;
282 grp = H5Gcreate(fid,
"Chombo_global", H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
287 int ncomp = num_output_real + num_output_int;
294 int real_comp_count = AMREX_SPACEDIM;
296 for (
int i = 0; i < NStructReal + pc.NumRealComps(); ++i ) {
297 if (write_real_comp[i]) {
299 snprintf(comp_name,
sizeof comp_name,
"real_component_%d", real_comp_count);
309 int int_comp_count = 2;
311 for (
int i = 0; i < NStructInt + pc.NumIntComps(); ++i ) {
312 if (write_int_comp[i]) {
313 snprintf(comp_name,
sizeof comp_name,
"int_component_%d", int_comp_count);
326 int finest_level = pc.finestLevel();
329 char level_name[128];
332 hid_t dcpl_id = H5Pcreate(H5P_DATASET_CREATE);
333 H5Pset_fill_time(dcpl_id, H5D_FILL_TIME_NEVER);
335 for (
int lev = 0; lev <= finest_level; ++lev) {
336 snprintf(level_name,
sizeof level_name,
"level_%d", lev);
338 grp = H5Gcreate(fid, level_name, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
340 std::cout <<
"H5Gcreate [" << level_name <<
"] failed!" << std::endl;
345 ngrids = pc.ParticleBoxArray(lev).size();
351 int mfs_size = 2 * AMREX_SPACEDIM;
352 hsize_t mfs_dim = (hsize_t)ngrids;
354 hid_t mfs_dset_space = H5Screate_simple(1, &mfs_dim, NULL);
355 hid_t mfs_dset= H5Dcreate(grp,
"boxes", comp_dtype, mfs_dset_space, H5P_DEFAULT, dcpl_id, H5P_DEFAULT);
357 Vector<int> vbox(ngrids * mfs_size);
358 for(
int j = 0; j < pc.ParticleBoxArray(lev).
size(); ++j) {
359 for(
int i = 0; i < AMREX_SPACEDIM; ++i) {
360 vbox[(j * mfs_size) + i] = pc.ParticleBoxArray(lev)[j].smallEnd(i);
361 vbox[(j * mfs_size) + i + AMREX_SPACEDIM] = pc.ParticleBoxArray(lev)[j].bigEnd(i);
365 status = H5Dwrite(mfs_dset, comp_dtype, H5S_ALL, H5S_ALL, H5P_DEFAULT, &(vbox[0]));
367 std::string msg(
"ParticleContainer::WriteHDF5ParticleDataSync(): unable to write boxes dataset");
371 H5Sclose(mfs_dset_space);
387 ParmParse
pp(
"particles");
388 pp.
query(
"particles_nfiles",nOutFiles);
389 if(nOutFiles == -1) nOutFiles =
NProcs;
391 pc.nOutFilesPrePost = nOutFiles;
393 fapl = H5Pcreate (H5P_FILE_ACCESS);
400 #ifdef AMREX_USE_HDF5_ASYNC
401 fid = H5Fopen_async(HDF5FileName.c_str(), H5F_ACC_RDWR, fapl, es_par_g);
403 fid = H5Fopen(HDF5FileName.c_str(), H5F_ACC_RDWR, fapl);
410 for (
int lev = 0; lev <= pc.finestLevel(); lev++)
412 snprintf(level_name,
sizeof level_name,
"level_%d", lev);
413 #ifdef AMREX_USE_HDF5_ASYNC
414 grp = H5Gopen_async(fid, level_name, H5P_DEFAULT, es_par_g);
416 grp = H5Gopen(fid, level_name, H5P_DEFAULT);
421 gotsome = (pc.nParticlesAtLevelPrePost[lev] > 0);
423 gotsome = (pc.NumberOfParticlesAtLevel(lev) > 0);
427 info.SetAlloc(
false);
428 MultiFab state(pc.ParticleBoxArray(lev),
429 pc.ParticleDistributionMap(lev),
434 Vector<int> which(state.size(),0);
435 Vector<int > count(state.size(),0);
436 Vector<Long> where(state.size(),0);
439 pc.filePrefixPrePost[lev] = HDF5FileName;
444 pc.WriteParticlesHDF5(lev, grp, which, count, where,
445 write_real_comp, write_int_comp,
447 particle_io_flags, is_checkpoint);
450 pc.whichPrePost[lev] = which;
451 pc.countPrePost[lev] = count;
452 pc.wherePrePost[lev] = where;
460 #ifdef AMREX_USE_HDF5_ASYNC
461 H5Gclose_async(grp, es_par_g);
467 H5Tclose(comp_dtype);
470 #ifdef AMREX_USE_HDF5_ASYNC
471 H5Fclose_async(fid, es_par_g);
#define BL_PROFILE(a)
Definition: AMReX_BLProfiler.H:551
#define BL_PROFILE_VAR_STOP(vname)
Definition: AMReX_BLProfiler.H:563
#define BL_PROFILE_VAR(fname, vname)
Definition: AMReX_BLProfiler.H:560
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
#define AMREX_ALWAYS_ASSERT(EX)
Definition: AMReX_BLassert.H:50
amrex::ParmParse pp
Input file parser instance for the given namespace.
Definition: AMReX_HypreIJIface.cpp:15
static int CreateWriteHDF5Attr(hid_t loc, const char *name, hsize_t n, void *data, hid_t dtype)
Definition: AMReX_WriteBinaryParticleDataHDF5.H:4
static int ReadHDF5Attr(hid_t loc, const char *name, void *data, hid_t dtype)
Definition: AMReX_WriteBinaryParticleDataHDF5.H:59
static int CreateWriteHDF5AttrString(hid_t loc, const char *name, const char *str)
Definition: AMReX_WriteBinaryParticleDataHDF5.H:28
static void SetHDF5fapl(hid_t fapl, MPI_Comm comm)
Definition: AMReX_WriteBinaryParticleDataHDF5.H:80
void WriteHDF5ParticleDataSync(PC const &pc, const std::string &dir, const std::string &name, const Vector< int > &write_real_comp, const Vector< int > &write_int_comp, const Vector< std::string > &real_comp_names, const Vector< std::string > &int_comp_names, const std::string &compression, F &&f, bool is_checkpoint)
Definition: AMReX_WriteBinaryParticleDataHDF5.H:115
int MPI_Comm
Definition: AMReX_ccse-mpi.H:47
int query(const char *name, bool &ref, int ival=FIRST) const
Same as querykth() but searches for the last occurrence of name.
Definition: AMReX_ParmParse.cpp:1309
This class provides the user with a few print options.
Definition: AMReX_Print.H:35
AMREX_GPU_HOST_DEVICE Long size(T const &b) noexcept
integer version
Definition: AMReX_GpuRange.H:26
void streamSynchronize() noexcept
Definition: AMReX_GpuDevice.H:237
int NProcs()
Process ID in MPI_COMM_WORLD.
Definition: AMReX_MPMD.cpp:122
MPI_Comm Communicator() noexcept
Definition: AMReX_ParallelDescriptor.H:210
void ReduceIntSum(int &)
Integer sum reduction.
Definition: AMReX_ParallelDescriptor.cpp:1252
void ReduceLongMax(Long &)
Long max reduction.
Definition: AMReX_ParallelDescriptor.cpp:1224
void ReduceLongSum(Long &)
Long sum reduction.
Definition: AMReX_ParallelDescriptor.cpp:1223
int IOProcessorNumber() noexcept
Definition: AMReX_ParallelDescriptor.H:266
void Barrier(const std::string &)
Definition: AMReX_ParallelDescriptor.cpp:1202
bool IOProcessor() noexcept
Is this CPU the I/O Processor? To get the rank number, call IOProcessorNumber()
Definition: AMReX_ParallelDescriptor.H:275
static int f(amrex::Real t, N_Vector y_data, N_Vector y_rhs, void *user_data)
Definition: AMReX_SundialsIntegrator.H:44
void FileOpenFailed(const std::string &file)
Output a message and abort when couldn't open the file.
Definition: AMReX_Utility.cpp:131
void CreateDirectoryFailed(const std::string &dir)
Output a message and abort when couldn't create the directory.
Definition: AMReX_Utility.cpp:123
bool UtilCreateDirectory(const std::string &path, mode_t mode, bool verbose=false)
Creates the specified directories. path may be either a full pathname or a relative pathname....
Definition: AMReX_Utility.cpp:116
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)
This shuts up the compiler about unused variables.
Definition: AMReX.H:111
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:225
void ExecOnFinalize(std::function< void()>)
We maintain a stack of functions that need to be called in Finalize(). The functions are called in LI...
Definition: AMReX.cpp:314
std::enable_if_t< RunOnGpu< typename Container< int, Allocator >::allocator_type >::value > fillFlags(Container< int, Allocator > &pflags, const PTile &ptile, F const &f)
Definition: AMReX_WriteBinaryParticleData.H:34
std::enable_if_t< RunOnGpu< typename Container< int, Allocator >::allocator_type >::value, amrex::Long > countFlags(const Vector< std::map< std::pair< int, int >, Container< int, Allocator >>> &particle_io_flags, const PC &pc)
Definition: AMReX_WriteBinaryParticleData.H:76