Skip to content
Snippets Groups Projects
Commit c3b8c5ca authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Also read/write star particles when using parallel HDF5

parent 9d330d2e
No related branches found
No related tags found
1 merge request!310Star particles and gparts links over MPI
......@@ -384,9 +384,10 @@ int main(int argc, char *argv[]) {
if (myrank == 0) clocks_gettime(&tic);
#if defined(WITH_MPI)
#if defined(HAVE_PARALLEL_HDF5)
read_ic_parallel(ICfileName, &us, dim, &parts, &gparts, &Ngas, &Ngpart,
&periodic, &flag_entropy_ICs, myrank, nr_nodes,
MPI_COMM_WORLD, MPI_INFO_NULL, dry_run);
read_ic_parallel(ICfileName, &us, dim, &parts, &gparts, &sparts, &Ngas, &Ngpart,
&Nspart, &periodic, &flag_entropy_ICs, with_hydro,
(with_external_gravity || with_self_gravity), with_stars,
myrank, nr_nodes, MPI_COMM_WORLD, MPI_INFO_NULL, dry_run);
#else
read_ic_serial(ICfileName, &us, dim, &parts, &gparts, &sparts, &Ngas, &Ngpart,
&Nspart, &periodic, &flag_entropy_ICs, with_hydro,
......
......@@ -46,6 +46,7 @@
#include "io_properties.h"
#include "kernel_hydro.h"
#include "part.h"
#include "stars_io.h"
#include "units.h"
/**
......@@ -373,9 +374,12 @@ void writeArray(struct engine* e, hid_t grp, char* fileName, FILE* xmfFile,
*/
void read_ic_parallel(char* fileName, const struct UnitSystem* internal_units,
double dim[3], struct part** parts, struct gpart** gparts,
size_t* Ngas, size_t* Ngparts, int* periodic,
int* flag_entropy, int mpi_rank, int mpi_size,
MPI_Comm comm, MPI_Info info, int dry_run) {
struct spart** sparts, size_t* Ngas, size_t* Ngparts,
size_t* Nstars, int* periodic, int* flag_entropy,
int with_hydro, int with_gravity, int with_stars,
int mpi_rank, int mpi_size, MPI_Comm comm, MPI_Info info,
int dry_run) {
hid_t h_file = 0, h_grp = 0;
/* GADGET has only cubic boxes (in cosmological mode) */
double boxSize[3] = {0.0, -1.0, -1.0};
......@@ -385,6 +389,7 @@ void read_ic_parallel(char* fileName, const struct UnitSystem* internal_units,
long long N_total[NUM_PARTICLE_TYPES] = {0};
long long offset[NUM_PARTICLE_TYPES] = {0};
int dimension = 3; /* Assume 3D if nothing is specified */
size_t Ndm;
/* Open file */
/* message("Opening file '%s' as IC.", fileName); */
......@@ -492,29 +497,38 @@ void read_ic_parallel(char* fileName, const struct UnitSystem* internal_units,
units_conversion_factor(ic_units, internal_units, UNIT_CONV_LENGTH);
/* Allocate memory to store SPH particles */
*Ngas = N[0];
if (posix_memalign((void*)parts, part_align, (*Ngas) * sizeof(struct part)) !=
0)
error("Error while allocating memory for particles");
bzero(*parts, *Ngas * sizeof(struct part));
/* Allocate memory to store all particles */
const size_t Ndm = N[1];
*Ngparts = N[1] + N[0];
if (posix_memalign((void*)gparts, gpart_align,
*Ngparts * sizeof(struct gpart)) != 0)
error(
"Error while allocating memory for gravity "
"particles");
bzero(*gparts, *Ngparts * sizeof(struct gpart));
/* message("Allocated %8.2f MB for particles.", *N *
* sizeof(struct part) /
if (with_hydro) {
*Ngas = N[0];
if (posix_memalign((void*)parts, part_align, (*Ngas) * sizeof(struct part)) !=
0)
error("Error while allocating memory for particles");
bzero(*parts, *Ngas * sizeof(struct part));
}
/* Allocate memory to store star particles */
if (with_stars) {
*Nstars = N[STAR];
if (posix_memalign((void*)sparts, spart_align,
*Nstars * sizeof(struct spart)) != 0)
error("Error while allocating memory for star particles");
bzero(*sparts, *Nstars * sizeof(struct spart));
}
/* Allocate memory to store gravity particles */
if (with_gravity) {
Ndm = N[1];
*Ngparts = (with_hydro ? N[GAS] : 0) + N[DM] + (with_stars ? N[STAR] : 0);
if (posix_memalign((void*)gparts, gpart_align,
*Ngparts * sizeof(struct gpart)) != 0)
error("Error while allocating memory for gravity particles");
bzero(*gparts, *Ngparts * sizeof(struct gpart));
}
/* message("Allocated %8.2f MB for particles.", *N * sizeof(struct part) /
* (1024.*1024.)); */
/* message("BoxSize = %lf", dim[0]); */
/* message("NumPart = [%zd, %zd] Total = %zd", *Ngas, Ndm,
* *Ngparts); */
/* message("NumPart = [%zd, %zd] Total = %zd", *Ngas, Ndm, *Ngparts); */
/* Loop over all particle types */
for (int ptype = 0; ptype < NUM_PARTICLE_TYPES; ptype++) {
......@@ -539,15 +553,26 @@ void read_ic_parallel(char* fileName, const struct UnitSystem* internal_units,
switch (ptype) {
case GAS:
Nparticles = *Ngas;
hydro_read_particles(*parts, list, &num_fields);
if (with_hydro) {
Nparticles = *Ngas;
hydro_read_particles(*parts, list, &num_fields);
}
break;
case DM:
Nparticles = Ndm;
darkmatter_read_particles(*gparts, list, &num_fields);
if (with_gravity) {
Nparticles = Ndm;
darkmatter_read_particles(*gparts, list, &num_fields);
}
break;
case STAR:
if (with_stars) {
Nparticles = *Nstars;
star_read_particles(*sparts, list, &num_fields);
}
break;
default:
message("Particle Type %d not yet supported. Particles ignored", ptype);
}
......@@ -563,10 +588,15 @@ void read_ic_parallel(char* fileName, const struct UnitSystem* internal_units,
}
/* Prepare the DM particles */
if (!dry_run) prepare_dm_gparts(*gparts, Ndm);
if (!dry_run && with_gravity) prepare_dm_gparts(*gparts, Ndm);
/* Now duplicate the hydro particle into gparts */
if (!dry_run) duplicate_hydro_gparts(*parts, *gparts, *Ngas, Ndm);
/* Duplicate the hydro particles into gparts */
if (!dry_run && with_gravity && with_hydro)
duplicate_hydro_gparts(*parts, *gparts, *Ngas, Ndm);
/* Duplicate the star particles into gparts */
if (!dry_run && with_gravity && with_stars)
duplicate_star_gparts(*sparts, *gparts, *Nstars, Ndm + *Ngas);
/* message("Done Reading particles..."); */
......@@ -609,17 +639,19 @@ void write_output_parallel(struct engine* e, const char* baseName,
hid_t h_file = 0, h_grp = 0;
const size_t Ngas = e->s->nr_parts;
const size_t Nstars = e->s->nr_sparts;
const size_t Ntot = e->s->nr_gparts;
int periodic = e->s->periodic;
int numFiles = 1;
struct part* parts = e->s->parts;
struct gpart* gparts = e->s->gparts;
struct gpart* dmparts = NULL;
struct spart* sparts = e->s->sparts;
static int outputCount = 0;
FILE* xmfFile = 0;
/* Number of unassociated gparts */
const size_t Ndm = Ntot > 0 ? Ntot - Ngas : 0;
const size_t Ndm = Ntot > 0 ? Ntot - (Ngas + Nstars) : 0;
/* File name */
char fileName[FILENAME_BUFFER_SIZE];
......@@ -642,7 +674,7 @@ void write_output_parallel(struct engine* e, const char* baseName,
/* Compute offset in the file and total number of
* particles */
size_t N[NUM_PARTICLE_TYPES] = {Ngas, Ndm, 0};
size_t N[NUM_PARTICLE_TYPES] = {Ngas, Ndm, 0, 0, Nstars, 0};
long long N_total[NUM_PARTICLE_TYPES] = {0};
long long offset[NUM_PARTICLE_TYPES] = {0};
MPI_Exscan(&N, &offset, NUM_PARTICLE_TYPES, MPI_LONG_LONG_INT, MPI_SUM, comm);
......@@ -816,11 +848,13 @@ void write_output_parallel(struct engine* e, const char* baseName,
/* Write DM particles */
Nparticles = Ndm;
darkmatter_write_particles(dmparts, list, &num_fields);
/* Free temporary array */
free(dmparts);
break;
case STAR:
Nparticles = Nstars;
star_write_particles(sparts, list, &num_fields);
break;
default:
error("Particle Type %d not yet supported. Aborting", ptype);
}
......@@ -830,9 +864,12 @@ void write_output_parallel(struct engine* e, const char* baseName,
writeArray(e, h_grp, fileName, xmfFile, partTypeGroupName, list[i],
Nparticles, N_total[ptype], mpi_rank, offset[ptype],
internal_units, snapshot_units);
/* Free temporary array */
free(dmparts);
if (dmparts) {
free(dmparts);
dmparts = 0;
}
/* Close particle group */
H5Gclose(h_grp);
......
......@@ -36,9 +36,11 @@
void read_ic_parallel(char* fileName, const struct UnitSystem* internal_units,
double dim[3], struct part** parts, struct gpart** gparts,
size_t* Ngas, size_t* Ngparts, int* periodic,
int* flag_entropy, int mpi_rank, int mpi_size,
MPI_Comm comm, MPI_Info info, int dry_run);
struct spart** sparts, size_t* Ngas, size_t* Ngparts,
size_t* Nsparts, int* periodic, int* flag_entropy,
int with_hydro, int with_gravity, int with_stars,
int mpi_rank, int mpi_size, MPI_Comm comm, MPI_Info info,
int dry_run);
void write_output_parallel(struct engine* e, const char* baseName,
const struct UnitSystem* internal_units,
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment