Skip to content
Snippets Groups Projects
Commit 29e3855a authored by Cristian Lalescu's avatar Cristian Lalescu
Browse files

Merge branch 'feature/efficient-particle-sampler' into develop

parents 920f623e c5f7042b
Branches
Tags
1 merge request!23WIP: Feature/use cmake
Pipeline #
...@@ -2,7 +2,6 @@ ...@@ -2,7 +2,6 @@
#include <cmath> #include <cmath>
#include "NSVEparticles.hpp" #include "NSVEparticles.hpp"
#include "scope_timer.hpp" #include "scope_timer.hpp"
#include "particles/particles_sampling.hpp"
template <typename rnumber> template <typename rnumber>
int NSVEparticles<rnumber>::initialize(void) int NSVEparticles<rnumber>::initialize(void)
...@@ -27,6 +26,13 @@ int NSVEparticles<rnumber>::initialize(void) ...@@ -27,6 +26,13 @@ int NSVEparticles<rnumber>::initialize(void)
"tracers0", "tracers0",
nparticles, nparticles,
tracers0_integration_steps); tracers0_integration_steps);
this->particles_sample_writer_mpi = new particles_output_sampling_hdf5<
long long int, double, 3, 3>(
MPI_COMM_WORLD,
this->ps->getGlobalNbParticles(),
(this->simname + "_particles.h5"),
"tracers0",
"position/0");
return EXIT_SUCCESS; return EXIT_SUCCESS;
} }
...@@ -60,6 +66,7 @@ int NSVEparticles<rnumber>::finalize(void) ...@@ -60,6 +66,7 @@ int NSVEparticles<rnumber>::finalize(void)
{ {
this->ps.release(); this->ps.release();
delete this->particles_output_writer_mpi; delete this->particles_output_writer_mpi;
delete this->particles_sample_writer_mpi;
this->NSVE<rnumber>::finalize(); this->NSVE<rnumber>::finalize();
return EXIT_SUCCESS; return EXIT_SUCCESS;
} }
...@@ -77,30 +84,50 @@ int NSVEparticles<rnumber>::do_stats() ...@@ -77,30 +84,50 @@ int NSVEparticles<rnumber>::do_stats()
if (!(this->iteration % this->niter_part == 0)) if (!(this->iteration % this->niter_part == 0))
return EXIT_SUCCESS; return EXIT_SUCCESS;
// allocate temporary data array
std::unique_ptr<double[]> pdata(new double[3*this->ps->getLocalNbParticles()]);
// copy position data
/// sample position /// sample position
sample_particles_system_position( std::copy(this->ps->getParticlesPositions(),
this->ps, this->ps->getParticlesPositions()+3*this->ps->getLocalNbParticles(),
(this->simname + "_particles.h5"), // filename pdata.get());
"tracers0", // hdf5 parent group this->particles_sample_writer_mpi->save_dataset(
"position" // dataset basename TODO "tracers0",
); "position",
this->ps->getParticlesPositions(),
&pdata,
this->ps->getParticlesIndexes(),
this->ps->getLocalNbParticles(),
this->ps->get_step_idx()-1);
/// sample velocity /// sample velocity
sample_from_particles_system(*this->tmp_vec_field, // field to save this->ps->sample_compute_field(*this->tmp_vec_field, pdata.get());
this->ps, this->particles_sample_writer_mpi->save_dataset(
(this->simname + "_particles.h5"), // filename "tracers0",
"tracers0", // hdf5 parent group "velocity",
"velocity" // dataset basename TODO this->ps->getParticlesPositions(),
); &pdata,
this->ps->getParticlesIndexes(),
this->ps->getLocalNbParticles(),
this->ps->get_step_idx()-1);
/// compute acceleration and sample it /// compute acceleration and sample it
this->fs->compute_Lagrangian_acceleration(this->tmp_vec_field); this->fs->compute_Lagrangian_acceleration(this->tmp_vec_field);
this->tmp_vec_field->ift(); this->tmp_vec_field->ift();
sample_from_particles_system(*this->tmp_vec_field, this->ps->sample_compute_field(*this->tmp_vec_field, pdata.get());
this->ps, this->particles_sample_writer_mpi->save_dataset(
(this->simname + "_particles.h5"), "tracers0",
"tracers0", "acceleration",
"acceleration"); this->ps->getParticlesPositions(),
&pdata,
this->ps->getParticlesIndexes(),
this->ps->getLocalNbParticles(),
this->ps->get_step_idx()-1);
// deallocate temporary data array
pdata.release();
return EXIT_SUCCESS; return EXIT_SUCCESS;
} }
......
...@@ -35,6 +35,7 @@ ...@@ -35,6 +35,7 @@
#include "full_code/NSVE.hpp" #include "full_code/NSVE.hpp"
#include "particles/particles_system_builder.hpp" #include "particles/particles_system_builder.hpp"
#include "particles/particles_output_hdf5.hpp" #include "particles/particles_output_hdf5.hpp"
#include "particles/particles_sampling.hpp"
/** \brief Navier-Stokes solver that includes simple Lagrangian tracers. /** \brief Navier-Stokes solver that includes simple Lagrangian tracers.
* *
...@@ -58,6 +59,7 @@ class NSVEparticles: public NSVE<rnumber> ...@@ -58,6 +59,7 @@ class NSVEparticles: public NSVE<rnumber>
/* other stuff */ /* other stuff */
std::unique_ptr<abstract_particles_system<long long int, double>> ps; std::unique_ptr<abstract_particles_system<long long int, double>> ps;
particles_output_hdf5<long long int, double,3,3> *particles_output_writer_mpi; particles_output_hdf5<long long int, double,3,3> *particles_output_writer_mpi;
particles_output_sampling_hdf5<long long int, double, 3, 3> *particles_sample_writer_mpi;
NSVEparticles( NSVEparticles(
......
...@@ -24,7 +24,7 @@ int field_test<rnumber>::read_parameters() ...@@ -24,7 +24,7 @@ int field_test<rnumber>::read_parameters()
this->test::read_parameters(); this->test::read_parameters();
// in case any parameters are needed, this is where they should be read // in case any parameters are needed, this is where they should be read
hid_t parameter_file; hid_t parameter_file;
hid_t dset, memtype, space; hid_t dset;
parameter_file = H5Fopen( parameter_file = H5Fopen(
(this->simname + std::string(".h5")).c_str(), (this->simname + std::string(".h5")).c_str(),
H5F_ACC_RDONLY, H5F_ACC_RDONLY,
......
...@@ -40,7 +40,7 @@ int filter_test<rnumber>::read_parameters() ...@@ -40,7 +40,7 @@ int filter_test<rnumber>::read_parameters()
{ {
this->test::read_parameters(); this->test::read_parameters();
hid_t parameter_file; hid_t parameter_file;
hid_t dset, memtype, space; hid_t dset;
parameter_file = H5Fopen( parameter_file = H5Fopen(
(this->simname + std::string(".h5")).c_str(), (this->simname + std::string(".h5")).c_str(),
H5F_ACC_RDONLY, H5F_ACC_RDONLY,
......
...@@ -22,9 +22,8 @@ int test::main_loop(void) ...@@ -22,9 +22,8 @@ int test::main_loop(void)
int test::read_parameters() int test::read_parameters()
{ {
hid_t parameter_file; hid_t parameter_file;
hid_t dset, memtype, space; hid_t dset;
char fname[256]; char fname[256];
char *string_data;
sprintf(fname, "%s.h5", this->simname.c_str()); sprintf(fname, "%s.h5", this->simname.c_str());
parameter_file = H5Fopen(fname, H5F_ACC_RDONLY, H5P_DEFAULT); parameter_file = H5Fopen(fname, H5F_ACC_RDONLY, H5P_DEFAULT);
dset = H5Dopen(parameter_file, "/parameters/dealias_type", H5P_DEFAULT); dset = H5Dopen(parameter_file, "/parameters/dealias_type", H5P_DEFAULT);
......
...@@ -41,6 +41,10 @@ class abstract_particles_output { ...@@ -41,6 +41,10 @@ class abstract_particles_output {
partsize_t particles_chunk_current_offset; partsize_t particles_chunk_current_offset;
protected: protected:
MPI_Comm& getCom(){
return mpi_com;
}
MPI_Comm& getComWriter(){ MPI_Comm& getComWriter(){
return mpi_com_writer; return mpi_com_writer;
} }
......
...@@ -9,10 +9,11 @@ template <class partsize_t, ...@@ -9,10 +9,11 @@ template <class partsize_t,
class real_number, class real_number,
int size_particle_positions, int size_particle_positions,
int size_particle_rhs> int size_particle_rhs>
class particles_output_sampling_hdf5 : public abstract_particles_output<partsize_t, class particles_output_sampling_hdf5 : public abstract_particles_output<
real_number, partsize_t,
size_particle_positions, real_number,
size_particle_rhs>{ size_particle_positions,
size_particle_rhs>{
using Parent = abstract_particles_output<partsize_t, using Parent = abstract_particles_output<partsize_t,
real_number, real_number,
size_particle_positions, size_particle_positions,
...@@ -20,7 +21,7 @@ class particles_output_sampling_hdf5 : public abstract_particles_output<partsize ...@@ -20,7 +21,7 @@ class particles_output_sampling_hdf5 : public abstract_particles_output<partsize
hid_t file_id, pgroup_id; hid_t file_id, pgroup_id;
const std::string dataset_name; std::string dataset_name;
const bool use_collective_io; const bool use_collective_io;
public: public:
...@@ -34,7 +35,6 @@ public: ...@@ -34,7 +35,6 @@ public:
int dataset_exists = -1; int dataset_exists = -1;
if(my_rank == 0){ if(my_rank == 0){
// Parallel HDF5 write
hid_t file_id = H5Fopen( hid_t file_id = H5Fopen(
in_filename.c_str(), in_filename.c_str(),
H5F_ACC_RDWR | H5F_ACC_DEBUG, H5F_ACC_RDWR | H5F_ACC_DEBUG,
...@@ -54,16 +54,18 @@ public: ...@@ -54,16 +54,18 @@ public:
return dataset_exists; return dataset_exists;
} }
particles_output_sampling_hdf5(MPI_Comm in_mpi_com, particles_output_sampling_hdf5(
const partsize_t inTotalNbParticles, MPI_Comm in_mpi_com,
const std::string& in_filename, const partsize_t inTotalNbParticles,
const std::string& in_groupname, const std::string& in_filename,
const std::string& in_dataset_name, const std::string& in_groupname,
const bool in_use_collective_io = false) const std::string& in_dataset_name,
const bool in_use_collective_io = false)
: Parent(in_mpi_com, inTotalNbParticles, 1), : Parent(in_mpi_com, inTotalNbParticles, 1),
dataset_name(in_dataset_name), dataset_name(in_dataset_name),
use_collective_io(in_use_collective_io){ use_collective_io(in_use_collective_io){
if(Parent::isInvolved()){ if(Parent::isInvolved()){
// prepare parallel MPI access property list
hid_t plist_id_par = H5Pcreate(H5P_FILE_ACCESS); hid_t plist_id_par = H5Pcreate(H5P_FILE_ACCESS);
assert(plist_id_par >= 0); assert(plist_id_par >= 0);
int retTest = H5Pset_fapl_mpio( int retTest = H5Pset_fapl_mpio(
...@@ -72,7 +74,7 @@ public: ...@@ -72,7 +74,7 @@ public:
MPI_INFO_NULL); MPI_INFO_NULL);
assert(retTest >= 0); assert(retTest >= 0);
// Parallel HDF5 write // open file for parallel HDF5 access
file_id = H5Fopen( file_id = H5Fopen(
in_filename.c_str(), in_filename.c_str(),
H5F_ACC_RDWR | H5F_ACC_DEBUG, H5F_ACC_RDWR | H5F_ACC_DEBUG,
...@@ -81,6 +83,7 @@ public: ...@@ -81,6 +83,7 @@ public:
retTest = H5Pclose(plist_id_par); retTest = H5Pclose(plist_id_par);
assert(retTest >= 0); assert(retTest >= 0);
// open group
pgroup_id = H5Gopen( pgroup_id = H5Gopen(
file_id, file_id,
in_groupname.c_str(), in_groupname.c_str(),
...@@ -91,13 +94,65 @@ public: ...@@ -91,13 +94,65 @@ public:
~particles_output_sampling_hdf5(){ ~particles_output_sampling_hdf5(){
if(Parent::isInvolved()){ if(Parent::isInvolved()){
// close group
int retTest = H5Gclose(pgroup_id); int retTest = H5Gclose(pgroup_id);
assert(retTest >= 0); assert(retTest >= 0);
// close file
retTest = H5Fclose(file_id); retTest = H5Fclose(file_id);
assert(retTest >= 0); assert(retTest >= 0);
} }
} }
int switch_to_group(
const std::string &in_groupname)
{
if(Parent::isInvolved()){
// close old group
int retTest = H5Gclose(pgroup_id);
assert(retTest >= 0);
// open new group
pgroup_id = H5Gopen(
file_id,
in_groupname.c_str(),
H5P_DEFAULT);
assert(pgroup_id >= 0);
}
return EXIT_SUCCESS;
}
int save_dataset(
const std::string& in_groupname,
const std::string& in_dataset_name,
const real_number input_particles_positions[],
const std::unique_ptr<real_number[]> input_particles_rhs[],
const partsize_t index_particles[],
const partsize_t nb_particles,
const int idx_time_step)
{
// update group
int retTest = this->switch_to_group(
in_groupname);
assert(retTest == EXIT_SUCCESS);
// update dataset name
dataset_name = in_dataset_name + "/" + std::to_string(idx_time_step);
int dataset_exists;
if (this->getMyRank() == 0)
dataset_exists = H5Lexists(
pgroup_id,
dataset_name.c_str(),
H5P_DEFAULT);
AssertMpi(MPI_Bcast(&dataset_exists, 1, MPI_INT, 0, this->getCom()));
if (dataset_exists == 0)
this->save(
input_particles_positions,
input_particles_rhs,
index_particles,
nb_particles,
idx_time_step);
return EXIT_SUCCESS;
}
void write( void write(
const int /*idx_time_step*/, const int /*idx_time_step*/,
const real_number* /*particles_positions*/, const real_number* /*particles_positions*/,
...@@ -108,18 +163,26 @@ public: ...@@ -108,18 +163,26 @@ public:
TIMEZONE("particles_output_hdf5::write"); TIMEZONE("particles_output_hdf5::write");
assert(particles_idx_offset < Parent::getTotalNbParticles() || (particles_idx_offset == Parent::getTotalNbParticles() && nb_particles == 0)); assert(particles_idx_offset < Parent::getTotalNbParticles() ||
(particles_idx_offset == Parent::getTotalNbParticles() &&
nb_particles == 0));
assert(particles_idx_offset+nb_particles <= Parent::getTotalNbParticles()); assert(particles_idx_offset+nb_particles <= Parent::getTotalNbParticles());
static_assert(std::is_same<real_number, double>::value || static_assert(std::is_same<real_number, double>::value ||
std::is_same<real_number, float>::value, std::is_same<real_number, float>::value,
"real_number must be double or float"); "real_number must be double or float");
const hid_t type_id = (sizeof(real_number) == 8 ? H5T_NATIVE_DOUBLE : H5T_NATIVE_FLOAT); const hid_t type_id = (sizeof(real_number) == 8 ?
H5T_NATIVE_DOUBLE :
H5T_NATIVE_FLOAT);
hid_t plist_id = H5Pcreate(H5P_DATASET_XFER); hid_t plist_id = H5Pcreate(H5P_DATASET_XFER);
assert(plist_id >= 0); assert(plist_id >= 0);
{ {
int rethdf = H5Pset_dxpl_mpio(plist_id, use_collective_io ? H5FD_MPIO_COLLECTIVE : H5FD_MPIO_INDEPENDENT); int rethdf = H5Pset_dxpl_mpio(
plist_id,
(use_collective_io ?
H5FD_MPIO_COLLECTIVE :
H5FD_MPIO_INDEPENDENT));
assert(rethdf >= 0); assert(rethdf >= 0);
} }
{ {
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment