particles_output_sampling_hdf5.hpp 11.2 KB
Newer Older
1
2
3
4
/******************************************************************************
*                                                                             *
*  Copyright 2019 Max Planck Institute for Dynamics and Self-Organization     *
*                                                                             *
Cristian Lalescu's avatar
Cristian Lalescu committed
5
*  This file is part of TurTLE.                                               *
6
*                                                                             *
Cristian Lalescu's avatar
Cristian Lalescu committed
7
*  TurTLE is free software: you can redistribute it and/or modify             *
8
9
10
11
*  it under the terms of the GNU General Public License as published          *
*  by the Free Software Foundation, either version 3 of the License,          *
*  or (at your option) any later version.                                     *
*                                                                             *
Cristian Lalescu's avatar
Cristian Lalescu committed
12
*  TurTLE is distributed in the hope that it will be useful,                  *
13
14
15
16
17
*  but WITHOUT ANY WARRANTY; without even the implied warranty of             *
*  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the              *
*  GNU General Public License for more details.                               *
*                                                                             *
*  You should have received a copy of the GNU General Public License          *
Cristian Lalescu's avatar
Cristian Lalescu committed
18
*  along with TurTLE.  If not, see <http://www.gnu.org/licenses/>             *
19
20
21
22
23
24
25
*                                                                             *
* Contact: Cristian.Lalescu@ds.mpg.de                                         *
*                                                                             *
******************************************************************************/



26
27
28
29
30
31
32
33
34
#ifndef PARTICLES_OUTPUT_SAMPLING_HDF5_HPP
#define PARTICLES_OUTPUT_SAMPLING_HDF5_HPP

#include "abstract_particles_output.hpp"

#include <hdf5.h>

template <class partsize_t,
          class real_number,
35
          int size_particle_positions>
36
37
38
class particles_output_sampling_hdf5 : public abstract_particles_output<
                                       partsize_t,
                                       real_number,
39
                                       size_particle_positions>{
40
41
    using Parent = abstract_particles_output<partsize_t,
                                             real_number,
42
                                             size_particle_positions>;
43

44
45
    hid_t file_id, pgroup_id;

46
    std::string dataset_name;
47
    std::vector<hsize_t> particle_file_layout;   // to hold the shape of initial condition array
48
49
50
    const bool use_collective_io;

public:
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
    static bool DatasetExistsCol(MPI_Comm in_mpi_com,
                                  const std::string& in_filename,
                                  const std::string& in_groupname,
                                 const std::string& in_dataset_name){
        int my_rank;
        AssertMpi(MPI_Comm_rank(in_mpi_com, &my_rank));

        int dataset_exists = -1;

        if(my_rank == 0){
            hid_t file_id = H5Fopen(
                    in_filename.c_str(),
                    H5F_ACC_RDWR | H5F_ACC_DEBUG,
                    H5P_DEFAULT);
            assert(file_id >= 0);

            dataset_exists = H5Lexists(
                    file_id,
                    (in_groupname + "/" + in_dataset_name).c_str(),
                    H5P_DEFAULT);

            int retTest = H5Fclose(file_id);
            assert(retTest >= 0);
        }

        AssertMpi(MPI_Bcast( &dataset_exists, 1, MPI_INT, 0, in_mpi_com ));
        return dataset_exists;
    }

80
81
82
83
84
85
86
    particles_output_sampling_hdf5(
            MPI_Comm in_mpi_com,
            const partsize_t inTotalNbParticles,
            const std::string& in_filename,
            const std::string& in_groupname,
            const std::string& in_dataset_name,
            const bool in_use_collective_io = false)
87
88
            : Parent(in_mpi_com, inTotalNbParticles, 1),
              dataset_name(in_dataset_name),
89
              use_collective_io(in_use_collective_io){
90
        if(Parent::isInvolved()){
91
            // prepare parallel MPI access property list
92
93
94
95
96
97
            hid_t plist_id_par = H5Pcreate(H5P_FILE_ACCESS);
            assert(plist_id_par >= 0);
            int retTest = H5Pset_fapl_mpio(
                    plist_id_par,
                    Parent::getComWriter(),
                    MPI_INFO_NULL);
98
            variable_used_only_in_assert(retTest);
99
100
            assert(retTest >= 0);

101
            // open file for parallel HDF5 access
102
103
104
105
106
107
108
109
            file_id = H5Fopen(
                    in_filename.c_str(),
                    H5F_ACC_RDWR | H5F_ACC_DEBUG,
                    plist_id_par);
            assert(file_id >= 0);
            retTest = H5Pclose(plist_id_par);
            assert(retTest >= 0);

110
            // open group
111
112
113
114
115
116
            pgroup_id = H5Gopen(
                    file_id,
                    in_groupname.c_str(),
                    H5P_DEFAULT);
            assert(pgroup_id >= 0);
        }
117
118
119
    }

    ~particles_output_sampling_hdf5(){
120
        if(Parent::isInvolved()){
121
            // close group
122
            int retTest = H5Gclose(pgroup_id);
123
            variable_used_only_in_assert(retTest);
124
            assert(retTest >= 0);
125
            // close file
126
127
128
            retTest = H5Fclose(file_id);
            assert(retTest >= 0);
        }
129
    }
130

131
132
133
134
135
136
    int switch_to_group(
            const std::string &in_groupname)
    {
        if(Parent::isInvolved()){
            // close old group
            int retTest = H5Gclose(pgroup_id);
137
            variable_used_only_in_assert(retTest);
138
139
140
141
142
143
144
145
146
147
148
149
            assert(retTest >= 0);

            // open new group
            pgroup_id = H5Gopen(
                    file_id,
                    in_groupname.c_str(),
                    H5P_DEFAULT);
            assert(pgroup_id >= 0);
        }
        return EXIT_SUCCESS;
    }

150
    template <int size_particle_rhs>
151
152
153
154
155
156
157
158
159
160
161
162
    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);
163
        variable_used_only_in_assert(retTest);
164
165
166
167
168
169
170
171
172
173
174
        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)
175
            this->template save<size_particle_rhs>(
176
177
178
179
180
181
182
183
                input_particles_positions,
                input_particles_rhs,
                index_particles,
                nb_particles,
                idx_time_step);
        return EXIT_SUCCESS;
    }

184
    void write(
185
            const int /*idx_time_step*/,
186
187
188
            const real_number* /*particles_positions*/,
            const std::unique_ptr<real_number[]>* particles_rhs,
            const partsize_t nb_particles,
189
190
            const partsize_t particles_idx_offset,
            const int size_particle_rhs) final{
191
192
        assert(Parent::isInvolved());

193
        TIMEZONE("particles_output_sampling_hdf5::write");
194

195
196
197
        assert(particles_idx_offset < Parent::getTotalNbParticles() ||
               (particles_idx_offset == Parent::getTotalNbParticles() &&
                nb_particles == 0));
198
199
200
201
202
        assert(particles_idx_offset+nb_particles <= Parent::getTotalNbParticles());

        static_assert(std::is_same<real_number, double>::value ||
                      std::is_same<real_number, float>::value,
                      "real_number must be double or float");
203
204
205
        const hid_t type_id = (sizeof(real_number) == 8 ?
                               H5T_NATIVE_DOUBLE :
                               H5T_NATIVE_FLOAT);
206
207
208
209

        hid_t plist_id = H5Pcreate(H5P_DATASET_XFER);
        assert(plist_id >= 0);
        {
210
211
212
213
214
            int rethdf = H5Pset_dxpl_mpio(
                    plist_id,
                    (use_collective_io ?
                     H5FD_MPIO_COLLECTIVE :
                     H5FD_MPIO_INDEPENDENT));
215
            variable_used_only_in_assert(rethdf);
216
217
218
219
            assert(rethdf >= 0);
        }
        {
            assert(size_particle_rhs >= 0);
220
            std::vector<hsize_t> datacount = std::vector<hsize_t>(this->particle_file_layout);
221
            datacount.push_back(size_particle_rhs);
222
            hid_t dataspace = H5Screate_simple(datacount.size(), &datacount.front(), NULL);
223
224
            assert(dataspace >= 0);

225
            hid_t dataset_id = H5Dcreate( pgroup_id,
226
227
228
229
230
231
232
233
234
                                          dataset_name.c_str(),
                                          type_id,
                                          dataspace,
                                          H5P_DEFAULT,
                                          H5P_DEFAULT,
                                          H5P_DEFAULT);
            assert(dataset_id >= 0);

            assert(particles_idx_offset >= 0);
235
            const hsize_t count[2] = {
236
237
                hsize_t(nb_particles),
                hsize_t(size_particle_rhs)};
238
            const hsize_t offset[2] = {
239
240
                hsize_t(particles_idx_offset),
                0};
241
            hid_t memspace = H5Screate_simple(2, count, NULL);
242
243
            assert(memspace >= 0);

Cristian Lalescu's avatar
Cristian Lalescu committed
244
            const hsize_t file_count[2] = {hsize_t(Parent::getTotalNbParticles()), hsize_t(size_particle_rhs)};
245
            hid_t filespace = H5Screate_simple(2, file_count, NULL);
246
247
248
249
250
251
252
253
            assert(filespace >= 0);
            int rethdf = H5Sselect_hyperslab(
                    filespace,
                    H5S_SELECT_SET,
                    offset,
                    NULL,
                    count,
                    NULL);
254
            variable_used_only_in_assert(rethdf);
255
256
257
258
259
260
261
262
263
            assert(rethdf >= 0);

            herr_t	status = H5Dwrite(
                    dataset_id,
                    type_id,
                    memspace,
                    filespace,
                    plist_id,
                    particles_rhs[0].get());
264
            variable_used_only_in_assert(status);
265
266
267
268
269
270
271
272
273
274
275
            assert(status >= 0);
            rethdf = H5Sclose(filespace);
            assert(rethdf >= 0);
            rethdf = H5Sclose(memspace);
            assert(rethdf >= 0);
            rethdf = H5Dclose(dataset_id);
            assert(rethdf >= 0);
        }

        {
            int rethdf = H5Pclose(plist_id);
276
            variable_used_only_in_assert(rethdf);
277
278
279
            assert(rethdf >= 0);
        }
    }
280

281
    int setParticleFileLayout(const std::vector<hsize_t> input_layout){
282
283
284
285
286
287
288
289
290
        this->particle_file_layout.resize(input_layout.size());
        for (unsigned int i=0; i<this->particle_file_layout.size(); i++)
            this->particle_file_layout[i] = input_layout[i];
        return EXIT_SUCCESS;
    }

    std::vector<hsize_t> getParticleFileLayout(void){
        return std::vector<hsize_t>(this->particle_file_layout);
    }
291
292
293
};

#endif