diff options
| -rw-r--r-- | src/lib/tng_io.c | 29 | ||||
| -rw-r--r-- | src/lib/tng_io.h | 35 | ||||
| -rw-r--r-- | src/tests/tng_parallel_read.c | 85 |
3 files changed, 128 insertions, 21 deletions
diff --git a/src/lib/tng_io.c b/src/lib/tng_io.c index 3da1ceb..fd20180 100644 --- a/src/lib/tng_io.c +++ b/src/lib/tng_io.c @@ -6228,11 +6228,12 @@ tng_function_status tng_trajectory_copy(tng_trajectory_t src, dest->medium_stride_length = src->medium_stride_length; dest->long_stride_length = src->long_stride_length; - dest->n_particle_data_blocks = src->n_particle_data_blocks; - dest->n_data_blocks = src->n_data_blocks; - - dest->non_tr_particle_data = src->non_tr_particle_data; - dest->non_tr_data = src->non_tr_data; + /* Currently the non trajectory data blocks are not copied since it + * can lead to problems when freeing memory in a parallel block. */ + dest->n_particle_data_blocks = 0; + dest->n_data_blocks = 0; + dest->non_tr_particle_data = 0; + dest->non_tr_data = 0; frame_set->first_frame = -1; frame_set->n_mapping_blocks = 0; @@ -6974,10 +6975,10 @@ tng_function_status tng_num_frame_sets_get(const tng_trajectory_t tng_data, } tng_function_status tng_current_frame_set_get - (const tng_trajectory_t tng_data, - tng_trajectory_frame_set_t frame_set) + (tng_trajectory_t tng_data, + tng_trajectory_frame_set_t *frame_set) { - frame_set = &tng_data->current_trajectory_frame_set; + *frame_set = &tng_data->current_trajectory_frame_set; return(TNG_SUCCESS); } @@ -7625,6 +7626,18 @@ tng_function_status tng_frame_set_prev_frame_set_file_pos_get return(TNG_SUCCESS); } +tng_function_status tng_frame_set_frame_range_get + (const tng_trajectory_t tng_data, + const tng_trajectory_frame_set_t frame_set, + int64_t *first_frame, + int64_t *last_frame) +{ + *first_frame = frame_set->first_frame; + *last_frame = *first_frame + frame_set->n_frames - 1; + + return(TNG_SUCCESS); +} + /** Translate from the particle numbering used in a frame set to the real * particle numbering - used in the molecule description. * @param frame_set is the frame_set containing the mappings to use. diff --git a/src/lib/tng_io.h b/src/lib/tng_io.h index a5af682..38c9d43 100644 --- a/src/lib/tng_io.h +++ b/src/lib/tng_io.h @@ -353,6 +353,10 @@ tng_function_status tng_trajectory_destroy_(tng_trajectory_t *tng_data_p) /** * @brief Copy a trajectory data container (dest is setup as well). + * @details This does not copy all data - only what is absolute necessary for + * parallel i/o. This can be used inside pragma omp for setting up a thread + * local copy of src. It can be freed (using tng_trajectory_destroy at the + * end of the parallel block. * @param src the original trajectory. * @param dest_p a pointer to memory to initialise as a trajectory. * @details Memory is allocated during initialisation. @@ -964,11 +968,11 @@ tng_function_status tng_num_frame_sets_get_ * @return TNG_SUCCESS (0) if successful. */ tng_function_status tng_current_frame_set_get - (const tng_trajectory_t tng_data, - tng_trajectory_frame_set_t frame_set); + (tng_trajectory_t tng_data, + tng_trajectory_frame_set_t *frame_set); tng_function_status tng_current_frame_set_get_ - (const tng_trajectory_t tng_data, - tng_trajectory_frame_set_t frame_set) + (tng_trajectory_t tng_data, + tng_trajectory_frame_set_t *frame_set) { return(tng_current_frame_set_get(tng_data, frame_set)); } @@ -1048,6 +1052,29 @@ tng_function_status tng_frame_set_prev_frame_set_file_pos_get_ } /** + * @brief Get the first and last frames of the frame set. + * @param tng_data is a trajectory data container. + * @param frame_set is the frame set of which to get the frame range. + * @param first_frame is set to the first frame of the frame set. + * @param last_frame is set to the last frame of the frame set. + * @return TNG_SUCCESS (0) if successful. + */ +tng_function_status tng_frame_set_frame_range_get + (const tng_trajectory_t tng_data, + const tng_trajectory_frame_set_t frame_set, + int64_t *first_frame, + int64_t *last_frame); +tng_function_status tng_frame_set_frame_range_get_ + (const tng_trajectory_t tng_data, + const tng_trajectory_frame_set_t frame_set, + int64_t *first_frame, + int64_t *last_frame) +{ + return(tng_frame_set_frame_range_get(tng_data, frame_set, first_frame, + last_frame)); +} + +/** * @brief Setup a molecule container. * @param tng_data is a trajectory data container. * @param molecule is the molecule to initialise. Memory must be preallocated. diff --git a/src/tests/tng_parallel_read.c b/src/tests/tng_parallel_read.c index 2aaa5f7..8078753 100644 --- a/src/tests/tng_parallel_read.c +++ b/src/tests/tng_parallel_read.c @@ -3,15 +3,24 @@ #include <tng_io.h> #include <omp.h> + +/* N.B. this code is for testing parallel reading of trajectory frame sets. The + * performance is not improved very much and is to a large extent limited by + * disk i/o. It can however be used as inspiration for writing parallel code + * using the TNG library. */ + int main(int argc, char **argv) { tng_trajectory_t traj, local_traj = 0; - union data_values ***positions = 0; // A 3-dimensional array to be populated + union data_values ***local_positions = 0; // A 3-dimensional array to be populated + union data_values **particle_pos = 0; int64_t n_particles, n_values_per_frame, n_frame_sets, n_frames; + int64_t n_frames_per_frame_set; tng_data_type data_type; - int i; - int64_t particle = 0; + int i, j; + int64_t particle = 0, local_first_frame, local_last_frame; char atom_name[64], res_name[64]; + tng_trajectory_frame_set_t frame_set = 0; if(argc <= 1) { @@ -30,6 +39,8 @@ int main(int argc, char **argv) } tng_input_file_set(traj, argv[1]); + tng_current_frame_set_get(traj, &frame_set); + // Read the file headers tng_file_headers_read(traj, TNG_USE_HASH); @@ -39,6 +50,15 @@ int main(int argc, char **argv) } tng_num_frame_sets_get(traj, &n_frame_sets); + tng_num_frames_per_frame_set_get(traj, &n_frames_per_frame_set); + + particle_pos = malloc(sizeof(union data_values **) * n_frame_sets * + n_frames_per_frame_set); + for(i = n_frame_sets * n_frames_per_frame_set; i--;) + { + /* Assume 3 values per frame even if it's not determined yet */ + particle_pos[i] = malloc(sizeof(union data_values) * 3); + } printf("%"PRId64" frame sets\n", n_frame_sets); @@ -58,10 +78,13 @@ int main(int argc, char **argv) #pragma omp parallel \ private (n_frames, n_particles, n_values_per_frame, \ - data_type, positions) \ -firstprivate (local_traj) + local_first_frame, local_last_frame, j) \ +firstprivate (local_traj, local_positions, frame_set)\ +shared(data_type, traj, n_frame_sets, particle_pos, particle, i)\ +default(none) { - positions = 0; + /* Each tng_trajectory_t keeps its own file pointers and i/o positions. + * Therefore there must be a copy for each thread. */ tng_trajectory_copy(traj, &local_traj); #pragma omp for for(i = 0; i < n_frame_sets; i++) @@ -70,22 +93,66 @@ firstprivate (local_traj) { printf("FAILED finding frame set %d!\n", i); } - if(tng_particle_data_get(local_traj, TNG_TRAJ_POSITIONS, &positions, + if(tng_particle_data_get(local_traj, TNG_TRAJ_POSITIONS, &local_positions, &n_frames, &n_particles, &n_values_per_frame, &data_type) != TNG_SUCCESS) { printf("FAILED getting particle data\n"); } + tng_current_frame_set_get(local_traj, &frame_set); + tng_frame_set_frame_range_get(local_traj, frame_set, &local_first_frame, &local_last_frame); +// printf("Frame %"PRId64"-%"PRId64":\n", local_first_frame, local_last_frame); // printf("%"PRId64" %"PRId64" %"PRId64"\n", n_frames, n_particles, n_values_per_frame); + for(j = 0; j < n_frames; j++) + { + particle_pos[local_first_frame + j][0] = local_positions[j][particle][0]; + particle_pos[local_first_frame + j][1] = local_positions[j][particle][1]; + particle_pos[local_first_frame + j][2] = local_positions[j][particle][2]; + } } + // Free memory - if(positions) + if(local_positions) { - tng_particle_data_values_free(local_traj, positions, n_frames, n_particles, + tng_particle_data_values_free(local_traj, local_positions, n_frames, n_particles, n_values_per_frame, data_type); } tng_trajectory_destroy(&local_traj); } + switch(data_type) + { + case TNG_INT_DATA: + for(i = 0; i < n_frame_sets * n_frames_per_frame_set; i++) + { + printf("\t%"PRId64"\t%"PRId64"\t%"PRId64"\n", particle_pos[i][0].i, + particle_pos[i][1].i, particle_pos[i][2].i); + } + break; + case TNG_FLOAT_DATA: + for(i = 0; i < n_frame_sets * n_frames_per_frame_set; i++) + { + printf("\t%f\t%f\t%f\n", particle_pos[i][0].f, + particle_pos[i][1].f, particle_pos[i][2].f); + } + break; + case TNG_DOUBLE_DATA: + for(i = 0; i < n_frame_sets * n_frames_per_frame_set; i++) + { + printf("\t%f\t%f\t%f\n", particle_pos[i][0].d, + particle_pos[i][1].d, particle_pos[i][2].d); + } + break; + default: + break; + } + + /* Free more memory */ + for(i = n_frame_sets * n_frames_per_frame_set; i--;) + { + free(particle_pos[i]); + } + free(particle_pos); + tng_trajectory_destroy(&traj); return(0); |