particles_output_hdf5.hpp 15 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
/******************************************************************************
*                                                                             *
*  Copyright 2019 Max Planck Institute for Dynamics and Self-Organization     *
*                                                                             *
*  This file is part of bfps.                                                 *
*                                                                             *
*  bfps is free software: you can redistribute it and/or modify               *
*  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.                                     *
*                                                                             *
*  bfps is distributed in the hope that it will be useful,                    *
*  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          *
*  along with bfps.  If not, see <http://www.gnu.org/licenses/>               *
*                                                                             *
* Contact: Cristian.Lalescu@ds.mpg.de                                         *
*                                                                             *
******************************************************************************/



26
27
28
29
30
31
#ifndef PARTICLES_OUTPUT_HDF5_HPP
#define PARTICLES_OUTPUT_HDF5_HPP

#include <memory>
#include <vector>
#include <hdf5.h>
32
#include <sys/stat.h>
33
34
35
36

#include "abstract_particles_output.hpp"
#include "scope_timer.hpp"

37
38
template <class partsize_t,
          class real_number,
39
          int size_particle_positions>
40
41
class particles_output_hdf5 : public abstract_particles_output<partsize_t,
                                                               real_number,
42
                                                               size_particle_positions>{
43
44
    using Parent = abstract_particles_output<partsize_t,
                                             real_number,
45
                                             size_particle_positions>;
46

47
    std::string particle_species_name;
48

49
    hid_t file_id;
50
    const partsize_t total_nb_particles;
51
    std::vector<hsize_t> particle_file_layout;   // to hold the shape of initial condition array
52

53
54
    hid_t dset_id_state;
    hid_t dset_id_rhs;
55

56
57
    bool use_collective_io;

58
public:
59
60
    particles_output_hdf5(MPI_Comm in_mpi_com,
                          const std::string ps_name,
61
                          const partsize_t inTotalNbParticles,
62
63
                          const int in_nb_rhs,
                          const bool in_use_collective_io = false)
64
65
            : abstract_particles_output<partsize_t,
                                        real_number,
66
                                        size_particle_positions>(
67
68
69
70
71
72
73
                                                in_mpi_com,
                                                inTotalNbParticles,
                                                in_nb_rhs),
              particle_species_name(ps_name),
              file_id(0),
              total_nb_particles(inTotalNbParticles),
              dset_id_state(0),
74
75
              dset_id_rhs(0),
              use_collective_io(in_use_collective_io){}
76

77
    int open_file(std::string filename){
78
79
        if(Parent::isInvolved()){
            TIMEZONE("particles_output_hdf5::open_file");
80

81
82
83
84
85
86
87
88
            this->require_checkpoint_groups(filename);

            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);
89
            variable_used_only_in_assert(retTest);
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
            assert(retTest >= 0);

            // Parallel HDF5 write
            file_id = H5Fopen(
                    filename.c_str(),
                    H5F_ACC_RDWR | H5F_ACC_DEBUG,
                    plist_id_par);
            // file_id = H5Fcreate(filename.c_str(), H5F_ACC_TRUNC | H5F_ACC_DEBUG/*H5F_ACC_EXCL*/, H5P_DEFAULT/*H5F_ACC_RDWR*/, plist_id_par);
            assert(file_id >= 0);
            H5Pclose(plist_id_par);

            dset_id_state = H5Gopen(
                    file_id,
                    (this->particle_species_name + std::string("/state")).c_str(),
                    H5P_DEFAULT);
            assert(dset_id_state >= 0);
            dset_id_rhs = H5Gopen(
                    file_id,
                    (this->particle_species_name + std::string("/rhs")).c_str(),
                    H5P_DEFAULT);
            assert(dset_id_rhs >= 0);
        }
112
        return EXIT_SUCCESS;
113
114
    }

115
    ~particles_output_hdf5(){}
116

117
118
119
120
121
122
    void update_particle_species_name(
            const std::string new_name)
    {
        this->particle_species_name.assign(new_name);
    }

123
    int close_file(void){
124
125
        if(Parent::isInvolved()){
            TIMEZONE("particles_output_hdf5::close_file");
126

127
            int rethdf = H5Gclose(dset_id_state);
128
            variable_used_only_in_assert(rethdf);
129
            assert(rethdf >= 0);
130

131
132
            rethdf = H5Gclose(dset_id_rhs);
            assert(rethdf >= 0);
133

134
135
136
            rethdf = H5Fclose(file_id);
            assert(rethdf >= 0);
        }
137
138
139
140
        return EXIT_SUCCESS;
    }

    void require_checkpoint_groups(std::string filename){
141
142
        if(Parent::isInvolved()){
            if (Parent::getMyRank() == 0)
143
            {
144
145
146
147
148
149
150
151
                bool file_exists = false;
                {
                    struct stat file_buffer;
                    file_exists = (stat(filename.c_str(), &file_buffer) == 0);
                }
                hid_t file_id;
                if (file_exists)
                    file_id = H5Fopen(
152
153
154
                        filename.c_str(),
                        H5F_ACC_RDWR | H5F_ACC_DEBUG,
                        H5P_DEFAULT);
155
156
157
158
159
160
                else
                    file_id = H5Fcreate(
                        filename.c_str(),
                        H5F_ACC_EXCL | H5F_ACC_DEBUG,
                        H5P_DEFAULT,
                        H5P_DEFAULT);
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
                assert(file_id >= 0);
                bool group_exists = H5Lexists(
                        file_id,
                        this->particle_species_name.c_str(),
                        H5P_DEFAULT);
                if (!group_exists)
                {
                    hid_t gg = H5Gcreate(
                        file_id,
                        this->particle_species_name.c_str(),
                        H5P_DEFAULT,
                        H5P_DEFAULT,
                        H5P_DEFAULT);
                    assert(gg >= 0);
                    H5Gclose(gg);
                }
                hid_t gg = H5Gopen(
                        file_id,
                        this->particle_species_name.c_str(),
                        H5P_DEFAULT);
181
                assert(gg >= 0);
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
                group_exists = H5Lexists(
                        gg,
                        "state",
                        H5P_DEFAULT);
                if (!group_exists)
                {
                    hid_t ggg = H5Gcreate(
                        gg,
                        "state",
                        H5P_DEFAULT,
                        H5P_DEFAULT,
                        H5P_DEFAULT);
                    assert(ggg >= 0);
                    H5Gclose(ggg);
                }
                group_exists = H5Lexists(
                        gg,
                        "rhs",
                        H5P_DEFAULT);
                if (!group_exists)
                {
                    hid_t ggg = H5Gcreate(
                        gg,
                        "rhs",
                        H5P_DEFAULT,
                        H5P_DEFAULT,
                        H5P_DEFAULT);
                    assert(ggg >= 0);
                    H5Gclose(ggg);
                }
212
                H5Gclose(gg);
213
                H5Fclose(file_id);
214
            }
215
            MPI_Barrier(Parent::getComWriter());
216
        }
217
218
    }

219
220
221
222
    void write(
            const int idx_time_step,
            const real_number* particles_positions,
            const std::unique_ptr<real_number[]>* particles_rhs,
223
            const partsize_t nb_particles,
224
225
            const partsize_t particles_idx_offset,
            const int size_particle_rhs) final{
226
227
        assert(Parent::isInvolved());

228
229
        TIMEZONE("particles_output_hdf5::write");

230
        assert(particles_idx_offset < Parent::getTotalNbParticles() || (particles_idx_offset == Parent::getTotalNbParticles() && nb_particles == 0));
231
232
        assert(particles_idx_offset+nb_particles <= Parent::getTotalNbParticles());

233
234
235
236
        static_assert(std::is_same<real_number, double>::value ||
                      std::is_same<real_number, float>::value,
                      "real_number must be double or float");
        const hid_t type_id = (sizeof(real_number) == 8 ? H5T_NATIVE_DOUBLE : H5T_NATIVE_FLOAT);
237

238
239
240
        hid_t plist_id = H5Pcreate(H5P_DATASET_XFER);
        assert(plist_id >= 0);
        {
241
            int rethdf = H5Pset_dxpl_mpio(plist_id, use_collective_io ? H5FD_MPIO_COLLECTIVE : H5FD_MPIO_INDEPENDENT);
242
            variable_used_only_in_assert(rethdf);
243
244
245
            assert(rethdf >= 0);
        }

246
        {
247
248
249
            std::vector<hsize_t> datacount = std::vector<hsize_t>(this->particle_file_layout);
            datacount.push_back(size_particle_positions);
            hid_t dataspace = H5Screate_simple(datacount.size(), &datacount.front(), NULL);
250
251
            assert(dataspace >= 0);

252
253
254
255
256
257
258
            hid_t dataset_id = H5Dcreate( dset_id_state,
                                          std::to_string(idx_time_step).c_str(),
                                          type_id,
                                          dataspace,
                                          H5P_DEFAULT,
                                          H5P_DEFAULT,
                                          H5P_DEFAULT);
259
260
            assert(dataset_id >= 0);

261
262
            assert(nb_particles >= 0);
            assert(particles_idx_offset >= 0);
263
264
265
            const hsize_t count[2] = {hsize_t(nb_particles), size_particle_positions};
            const hsize_t offset[2] = {hsize_t(particles_idx_offset), 0};
            hid_t memspace = H5Screate_simple(2, count, NULL);
266
267
            assert(memspace >= 0);

268
269
270
271
272
273
            assert(total_nb_particles >= 0);
            assert(size_particle_positions >= 0);
            const hsize_t file_count[2] = {hsize_t(total_nb_particles), size_particle_positions};
            hid_t filespace = H5Screate_simple(2, file_count, NULL);
            assert(filespace >= 0);

274
275
276
277
278
279
280
            int rethdf = H5Sselect_hyperslab(
                    filespace,
                    H5S_SELECT_SET,
                    offset,
                    NULL,
                    count,
                    NULL);
281
            variable_used_only_in_assert(rethdf);
282
283
            assert(rethdf >= 0);

284
285
286
287
288
289
290
            herr_t	status = H5Dwrite(
                    dataset_id,
                    type_id,
                    memspace,
                    filespace,
                    plist_id,
                    particles_positions);
291
            variable_used_only_in_assert(status);
292
293
294
295
296
297
298
299
300
            assert(status >= 0);
            rethdf = H5Sclose(memspace);
            assert(rethdf >= 0);
            rethdf = H5Dclose(dataset_id);
            assert(rethdf >= 0);
            rethdf = H5Sclose(filespace);
            assert(rethdf >= 0);
        }
        {
301
            assert(size_particle_rhs >= 0);
302
303
            std::vector<hsize_t> datacount = std::vector<hsize_t>(this->particle_file_layout);
            datacount.insert(datacount.begin(), hsize_t(Parent::getNbRhs()));
304
            datacount.push_back(size_particle_rhs);
305
            hid_t dataspace = H5Screate_simple(datacount.size(), &datacount.front(), NULL);
306
307
            assert(dataspace >= 0);

308
309
310
311
312
313
314
            hid_t dataset_id = H5Dcreate( dset_id_rhs,
                                          std::to_string(idx_time_step).c_str(),
                                          type_id,
                                          dataspace,
                                          H5P_DEFAULT,
                                          H5P_DEFAULT,
                                          H5P_DEFAULT);
315
316
            assert(dataset_id >= 0);

317
            assert(particles_idx_offset >= 0);
Berenger Bramas's avatar
Berenger Bramas committed
318
            for(int idx_rhs = 0 ; idx_rhs < Parent::getNbRhs() ; ++idx_rhs){
319
320
321
322
323
324
325
326
                const hsize_t count[3] = {
                    1,
                    hsize_t(nb_particles),
                    hsize_t(size_particle_rhs)};
                const hsize_t offset[3] = {
                    hsize_t(idx_rhs),
                    hsize_t(particles_idx_offset),
                    0};
Berenger Bramas's avatar
Berenger Bramas committed
327
328
329
                hid_t memspace = H5Screate_simple(3, count, NULL);
                assert(memspace >= 0);

330
331
332
333
                assert(total_nb_particles >= 0);
                assert(size_particle_positions >= 0);
                const hsize_t file_count[3] = {hsize_t(Parent::getNbRhs()), hsize_t(total_nb_particles), size_particle_positions};
                hid_t filespace = H5Screate_simple(3, file_count, NULL);
Berenger Bramas's avatar
Berenger Bramas committed
334
                assert(filespace >= 0);
335

336
337
338
339
340
341
342
                int rethdf = H5Sselect_hyperslab(
                        filespace,
                        H5S_SELECT_SET,
                        offset,
                        NULL,
                        count,
                        NULL);
343
                variable_used_only_in_assert(rethdf);
Berenger Bramas's avatar
Berenger Bramas committed
344
345
                assert(rethdf >= 0);

346
347
348
349
350
351
352
                herr_t	status = H5Dwrite(
                        dataset_id,
                        type_id,
                        memspace,
                        filespace,
                        plist_id,
                        particles_rhs[idx_rhs].get());
353
                variable_used_only_in_assert(status);
Berenger Bramas's avatar
Berenger Bramas committed
354
355
356
357
358
359
360
                assert(status >= 0);
                rethdf = H5Sclose(filespace);
                assert(rethdf >= 0);
                rethdf = H5Sclose(memspace);
                assert(rethdf >= 0);
            }
            int rethdf = H5Dclose(dataset_id);
361
            variable_used_only_in_assert(rethdf);
362
363
            assert(rethdf >= 0);
        }
364
365
366

        {
            int rethdf = H5Pclose(plist_id);
367
            variable_used_only_in_assert(rethdf);
368
369
            assert(rethdf >= 0);
        }
370
    }
371
372
373
374
375
376
377
378
379

    int setParticleFileLayout(std::vector<hsize_t> input_layout){
        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){
380
        return std::vector<hsize_t>(this->particle_file_layout);
381
    }
382
383
};

384
385
#endif//PARTICLES_OUTPUT_HDF5_HPP