From 37e526974d8af3304afbfcb5d408d5365bf4489b Mon Sep 17 00:00:00 2001 From: Benjamin Brands Date: Tue, 6 Feb 2018 07:59:04 +0100 Subject: [PATCH] enable chunking for save/load of ScaLAPACKMatrix --- include/deal.II/lac/scalapack.h | 16 +++- source/lac/scalapack.cc | 58 +++++++++-- tests/scalapack/scalapack_10_b.cc | 96 +++++++++++++++++++ ...calapack_10_b.with_hdf5=on.mpirun=1.output | 72 ++++++++++++++ ...alapack_10_b.with_hdf5=on.mpirun=11.output | 72 ++++++++++++++ ...calapack_10_b.with_hdf5=on.mpirun=4.output | 72 ++++++++++++++ ...calapack_10_b.with_hdf5=on.mpirun=5.output | 72 ++++++++++++++ 7 files changed, 446 insertions(+), 12 deletions(-) create mode 100644 tests/scalapack/scalapack_10_b.cc create mode 100644 tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=1.output create mode 100644 tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=11.output create mode 100644 tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=4.output create mode 100644 tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=5.output diff --git a/include/deal.II/lac/scalapack.h b/include/deal.II/lac/scalapack.h index 42468c81fc..9686e91918 100644 --- a/include/deal.II/lac/scalapack.h +++ b/include/deal.II/lac/scalapack.h @@ -191,6 +191,7 @@ public: /** * Stores the distributed matrix in @p filename using HDF5. + * * In case that deal.II was built without HDF5 * a call to this function will cause an exception to be thrown. * @@ -199,8 +200,15 @@ public: * internally the distributed matrix is copied to one process, which * does the output. Therefore, the matrix has to fit into the memory * of one process. + * + * To tweak the I/O performance, especially for parallel I/O, the user may define the optional parameter @p chunk_size. + * All MPI processes need to call the function with the same value. + * The matrix is written in chunks to the file, therefore the properties of the system define the optimal chunk size. + * Internally, HDF5 splits the matrix into chunk_size.first x chunk_size.second sized blocks, + * with chunk_size.first being the number of rows of a chunk and chunk_size.second the number of columns. */ - void save(const char *filename) const; + void save(const char *filename, + const std::pair &chunk_size=std::make_pair(numbers::invalid_unsigned_int,numbers::invalid_unsigned_int)) const; /** * Loads the distributed matrix from file @p filename using HDF5. @@ -417,7 +425,8 @@ private: * Stores the distributed matrix in @p filename * using serial routines */ - void save_serial(const char *filename) const; + void save_serial(const char *filename, + const std::pair &chunk_size) const; /* * Loads the distributed matrix from file @p filename @@ -429,7 +438,8 @@ private: * Stores the distributed matrix in @p filename * using parallel routines */ - void save_parallel(const char *filename) const; + void save_parallel(const char *filename, + const std::pair &chunk_size) const; /* * Loads the distributed matrix from file @p filename diff --git a/source/lac/scalapack.cc b/source/lac/scalapack.cc index 55d592b31a..f87838f1fa 100644 --- a/source/lac/scalapack.cc +++ b/source/lac/scalapack.cc @@ -914,19 +914,34 @@ NumberType ScaLAPACKMatrix::norm(const char type) const template -void ScaLAPACKMatrix::save(const char *filename) const +void ScaLAPACKMatrix::save(const char *filename, + const std::pair &chunk_size) const { #ifndef DEAL_II_WITH_HDF5 (void)filename; + (void)chunk_size; AssertThrow(false, ExcMessage ("HDF5 support is disabled.")); #else + + Assert((chunk_size.first <= (unsigned int)n_rows) && (chunk_size.first>0),ExcIndexRange(chunk_size.first,1,n_rows+1)); + Assert((chunk_size.second <= (unsigned int)n_columns) && (chunk_size.second>0),ExcIndexRange(chunk_size.second,1,n_columns+1)); + + std::pair chunks_size_ = chunk_size; + + if (chunks_size_.first==numbers::invalid_unsigned_int || chunks_size_.second==numbers::invalid_unsigned_int) + { + // default: store the matrix in chunks of columns + chunks_size_.first = n_rows; + chunks_size_.second = 1; + } + # ifdef H5_HAVE_PARALLEL //implementation for configurations equipped with a parallel file system - save_parallel(filename); + save_parallel(filename,chunks_size_); # else //implementation for configurations with no parallel file system - save_serial(filename); + save_serial(filename,chunks_size_); # endif #endif @@ -935,10 +950,12 @@ void ScaLAPACKMatrix::save(const char *filename) const template -void ScaLAPACKMatrix::save_serial(const char *filename) const +void ScaLAPACKMatrix::save_serial(const char *filename, + const std::pair &chunk_size) const { # ifndef DEAL_II_WITH_HDF5 (void)filename; + (void)chunk_size; Assert(false,ExcInternalError()); # else @@ -964,6 +981,15 @@ void ScaLAPACKMatrix::save_serial(const char *filename) const // create a new file using default properties hid_t file_id = H5Fcreate(filename, H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT); + // modify dataset creation properties, i.e. enable chunking + hsize_t chunk_dims[2]; + //revert order of rows and columns as ScaLAPACK uses column-major ordering + chunk_dims[0] = chunk_size.second; + chunk_dims[1] = chunk_size.first; + hid_t property = H5Pcreate (H5P_DATASET_CREATE); + status = H5Pset_chunk (property, 2, chunk_dims); + AssertThrow(status >= 0, ExcIO()); + // create the data space for the dataset hsize_t dims[2]; //change order of rows and columns as ScaLAPACKMatrix uses column major ordering @@ -971,11 +997,11 @@ void ScaLAPACKMatrix::save_serial(const char *filename) const dims[1] = n_rows; hid_t dataspace_id = H5Screate_simple(2, dims, nullptr); - // create the dataset + // create the dataset within the file using chunk creation properties hid_t type_id = hdf5_type_id(&tmp.values[0]); hid_t dataset_id = H5Dcreate2(file_id, "/matrix", type_id, dataspace_id, - H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); + H5P_DEFAULT, property, H5P_DEFAULT); // write the dataset status = H5Dwrite(dataset_id, type_id, @@ -991,6 +1017,10 @@ void ScaLAPACKMatrix::save_serial(const char *filename) const status = H5Sclose(dataspace_id); AssertThrow(status >= 0, ExcIO()); + // release the creation property + status = H5Pclose (property); + AssertThrow(status >= 0, ExcIO()); + // close the file. status = H5Fclose(file_id); AssertThrow(status >= 0, ExcIO()); @@ -1001,10 +1031,12 @@ void ScaLAPACKMatrix::save_serial(const char *filename) const template -void ScaLAPACKMatrix::save_parallel(const char *filename) const +void ScaLAPACKMatrix::save_parallel(const char *filename, + const std::pair &chunk_size) const { # ifndef DEAL_II_WITH_HDF5 (void)filename; + (void)chunk_size; Assert(false,ExcInternalError()); # else @@ -1048,12 +1080,20 @@ void ScaLAPACKMatrix::save_parallel(const char *filename) const hid_t filespace = H5Screate_simple(2, dims, nullptr); - // create the dataset with default properties and close filespace + // create the chunked dataset with default properties and close filespace + hsize_t chunk_dims[2]; + //revert order of rows and columns as ScaLAPACK uses column-major ordering + chunk_dims[0] = chunk_size.second; + chunk_dims[1] = chunk_size.first; + plist_id = H5Pcreate(H5P_DATASET_CREATE); + H5Pset_chunk(plist_id, 2, chunk_dims); hid_t type_id = hdf5_type_id(data); hid_t dset_id = H5Dcreate2(file_id, "/matrix", type_id, - filespace, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT); + filespace, H5P_DEFAULT, plist_id, H5P_DEFAULT); status = H5Sclose(filespace); AssertThrow(status >= 0, ExcIO()); + status = H5Pclose(plist_id); + AssertThrow(status >= 0, ExcIO()); // gather the number of local rows and columns from all processes std::vector proc_n_local_rows(n_mpi_processes), proc_n_local_columns(n_mpi_processes); diff --git a/tests/scalapack/scalapack_10_b.cc b/tests/scalapack/scalapack_10_b.cc new file mode 100644 index 0000000000..3e09b515bb --- /dev/null +++ b/tests/scalapack/scalapack_10_b.cc @@ -0,0 +1,96 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#include "../tests.h" +#include "../lapack/create_matrix.h" + +// test serial saving and loading of distributed ScaLAPACKMatrices with prescribed chunk sizes + +#include +#include +#include +#include +#include + +#include + +#include +#include +#include + + +template +void test(const std::pair &size, const unsigned int block_size, const std::pair &chunk_size) +{ + const std::string filename ("scalapck_10_b_test.h5"); + + MPI_Comm mpi_communicator(MPI_COMM_WORLD); + const unsigned int this_mpi_process(Utilities::MPI::this_mpi_process(mpi_communicator)); + ConditionalOStream pcout (std::cout, (this_mpi_process ==0)); + + FullMatrix full(size.first,size.second); + create_random(full); + + //create 2d process grid + std::shared_ptr grid = std::make_shared(mpi_communicator,size.first, + size.second,block_size,block_size); + + ScaLAPACKMatrix scalapack_matrix(size.first,size.second,grid,block_size,block_size); + ScaLAPACKMatrix scalapack_matrix_copy(size.first,size.second,grid,block_size,block_size); + + scalapack_matrix = full; + scalapack_matrix.save(filename.c_str(),chunk_size); + scalapack_matrix_copy.load(filename.c_str()); + + FullMatrix copy(size.first,size.second); + scalapack_matrix_copy.copy_to(copy); + copy.add(-1,full); + + pcout << size.first << "x" << size.second << " & " + << block_size << " & " + << chunk_size.first << "x" << chunk_size.second << std::endl; + AssertThrow(copy.frobenius_norm()<1e-12,ExcInternalError()); + std::remove(filename.c_str()); +} + + + +int main (int argc,char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, numbers::invalid_unsigned_int); + + std::vector> sizes; + sizes.push_back(std::make_pair(100,75)); + sizes.push_back(std::make_pair(200,225)); + sizes.push_back(std::make_pair(300,250)); + + const std::vector block_sizes = {{1,16,32}}; + + std::vector> chunk_sizes; + chunk_sizes.push_back(std::make_pair(1,1)); + chunk_sizes.push_back(std::make_pair(10,10)); + chunk_sizes.push_back(std::make_pair(50,50)); + chunk_sizes.push_back(std::make_pair(100,75)); + + for (unsigned int i=0; i(sizes[i],block_sizes[j],chunk_sizes[k]); + + for (unsigned int i=0; i(sizes[i],block_sizes[j],chunk_sizes[k]); +} diff --git a/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=1.output b/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=1.output new file mode 100644 index 0000000000..ec9538e8ac --- /dev/null +++ b/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=1.output @@ -0,0 +1,72 @@ +100x75 & 1 & 1x1 +100x75 & 1 & 10x10 +100x75 & 1 & 50x50 +100x75 & 1 & 100x75 +100x75 & 16 & 1x1 +100x75 & 16 & 10x10 +100x75 & 16 & 50x50 +100x75 & 16 & 100x75 +100x75 & 32 & 1x1 +100x75 & 32 & 10x10 +100x75 & 32 & 50x50 +100x75 & 32 & 100x75 +200x225 & 1 & 1x1 +200x225 & 1 & 10x10 +200x225 & 1 & 50x50 +200x225 & 1 & 100x75 +200x225 & 16 & 1x1 +200x225 & 16 & 10x10 +200x225 & 16 & 50x50 +200x225 & 16 & 100x75 +200x225 & 32 & 1x1 +200x225 & 32 & 10x10 +200x225 & 32 & 50x50 +200x225 & 32 & 100x75 +300x250 & 1 & 1x1 +300x250 & 1 & 10x10 +300x250 & 1 & 50x50 +300x250 & 1 & 100x75 +300x250 & 16 & 1x1 +300x250 & 16 & 10x10 +300x250 & 16 & 50x50 +300x250 & 16 & 100x75 +300x250 & 32 & 1x1 +300x250 & 32 & 10x10 +300x250 & 32 & 50x50 +300x250 & 32 & 100x75 +100x75 & 1 & 1x1 +100x75 & 1 & 10x10 +100x75 & 1 & 50x50 +100x75 & 1 & 100x75 +100x75 & 16 & 1x1 +100x75 & 16 & 10x10 +100x75 & 16 & 50x50 +100x75 & 16 & 100x75 +100x75 & 32 & 1x1 +100x75 & 32 & 10x10 +100x75 & 32 & 50x50 +100x75 & 32 & 100x75 +200x225 & 1 & 1x1 +200x225 & 1 & 10x10 +200x225 & 1 & 50x50 +200x225 & 1 & 100x75 +200x225 & 16 & 1x1 +200x225 & 16 & 10x10 +200x225 & 16 & 50x50 +200x225 & 16 & 100x75 +200x225 & 32 & 1x1 +200x225 & 32 & 10x10 +200x225 & 32 & 50x50 +200x225 & 32 & 100x75 +300x250 & 1 & 1x1 +300x250 & 1 & 10x10 +300x250 & 1 & 50x50 +300x250 & 1 & 100x75 +300x250 & 16 & 1x1 +300x250 & 16 & 10x10 +300x250 & 16 & 50x50 +300x250 & 16 & 100x75 +300x250 & 32 & 1x1 +300x250 & 32 & 10x10 +300x250 & 32 & 50x50 +300x250 & 32 & 100x75 diff --git a/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=11.output b/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=11.output new file mode 100644 index 0000000000..ec9538e8ac --- /dev/null +++ b/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=11.output @@ -0,0 +1,72 @@ +100x75 & 1 & 1x1 +100x75 & 1 & 10x10 +100x75 & 1 & 50x50 +100x75 & 1 & 100x75 +100x75 & 16 & 1x1 +100x75 & 16 & 10x10 +100x75 & 16 & 50x50 +100x75 & 16 & 100x75 +100x75 & 32 & 1x1 +100x75 & 32 & 10x10 +100x75 & 32 & 50x50 +100x75 & 32 & 100x75 +200x225 & 1 & 1x1 +200x225 & 1 & 10x10 +200x225 & 1 & 50x50 +200x225 & 1 & 100x75 +200x225 & 16 & 1x1 +200x225 & 16 & 10x10 +200x225 & 16 & 50x50 +200x225 & 16 & 100x75 +200x225 & 32 & 1x1 +200x225 & 32 & 10x10 +200x225 & 32 & 50x50 +200x225 & 32 & 100x75 +300x250 & 1 & 1x1 +300x250 & 1 & 10x10 +300x250 & 1 & 50x50 +300x250 & 1 & 100x75 +300x250 & 16 & 1x1 +300x250 & 16 & 10x10 +300x250 & 16 & 50x50 +300x250 & 16 & 100x75 +300x250 & 32 & 1x1 +300x250 & 32 & 10x10 +300x250 & 32 & 50x50 +300x250 & 32 & 100x75 +100x75 & 1 & 1x1 +100x75 & 1 & 10x10 +100x75 & 1 & 50x50 +100x75 & 1 & 100x75 +100x75 & 16 & 1x1 +100x75 & 16 & 10x10 +100x75 & 16 & 50x50 +100x75 & 16 & 100x75 +100x75 & 32 & 1x1 +100x75 & 32 & 10x10 +100x75 & 32 & 50x50 +100x75 & 32 & 100x75 +200x225 & 1 & 1x1 +200x225 & 1 & 10x10 +200x225 & 1 & 50x50 +200x225 & 1 & 100x75 +200x225 & 16 & 1x1 +200x225 & 16 & 10x10 +200x225 & 16 & 50x50 +200x225 & 16 & 100x75 +200x225 & 32 & 1x1 +200x225 & 32 & 10x10 +200x225 & 32 & 50x50 +200x225 & 32 & 100x75 +300x250 & 1 & 1x1 +300x250 & 1 & 10x10 +300x250 & 1 & 50x50 +300x250 & 1 & 100x75 +300x250 & 16 & 1x1 +300x250 & 16 & 10x10 +300x250 & 16 & 50x50 +300x250 & 16 & 100x75 +300x250 & 32 & 1x1 +300x250 & 32 & 10x10 +300x250 & 32 & 50x50 +300x250 & 32 & 100x75 diff --git a/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=4.output b/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=4.output new file mode 100644 index 0000000000..ec9538e8ac --- /dev/null +++ b/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=4.output @@ -0,0 +1,72 @@ +100x75 & 1 & 1x1 +100x75 & 1 & 10x10 +100x75 & 1 & 50x50 +100x75 & 1 & 100x75 +100x75 & 16 & 1x1 +100x75 & 16 & 10x10 +100x75 & 16 & 50x50 +100x75 & 16 & 100x75 +100x75 & 32 & 1x1 +100x75 & 32 & 10x10 +100x75 & 32 & 50x50 +100x75 & 32 & 100x75 +200x225 & 1 & 1x1 +200x225 & 1 & 10x10 +200x225 & 1 & 50x50 +200x225 & 1 & 100x75 +200x225 & 16 & 1x1 +200x225 & 16 & 10x10 +200x225 & 16 & 50x50 +200x225 & 16 & 100x75 +200x225 & 32 & 1x1 +200x225 & 32 & 10x10 +200x225 & 32 & 50x50 +200x225 & 32 & 100x75 +300x250 & 1 & 1x1 +300x250 & 1 & 10x10 +300x250 & 1 & 50x50 +300x250 & 1 & 100x75 +300x250 & 16 & 1x1 +300x250 & 16 & 10x10 +300x250 & 16 & 50x50 +300x250 & 16 & 100x75 +300x250 & 32 & 1x1 +300x250 & 32 & 10x10 +300x250 & 32 & 50x50 +300x250 & 32 & 100x75 +100x75 & 1 & 1x1 +100x75 & 1 & 10x10 +100x75 & 1 & 50x50 +100x75 & 1 & 100x75 +100x75 & 16 & 1x1 +100x75 & 16 & 10x10 +100x75 & 16 & 50x50 +100x75 & 16 & 100x75 +100x75 & 32 & 1x1 +100x75 & 32 & 10x10 +100x75 & 32 & 50x50 +100x75 & 32 & 100x75 +200x225 & 1 & 1x1 +200x225 & 1 & 10x10 +200x225 & 1 & 50x50 +200x225 & 1 & 100x75 +200x225 & 16 & 1x1 +200x225 & 16 & 10x10 +200x225 & 16 & 50x50 +200x225 & 16 & 100x75 +200x225 & 32 & 1x1 +200x225 & 32 & 10x10 +200x225 & 32 & 50x50 +200x225 & 32 & 100x75 +300x250 & 1 & 1x1 +300x250 & 1 & 10x10 +300x250 & 1 & 50x50 +300x250 & 1 & 100x75 +300x250 & 16 & 1x1 +300x250 & 16 & 10x10 +300x250 & 16 & 50x50 +300x250 & 16 & 100x75 +300x250 & 32 & 1x1 +300x250 & 32 & 10x10 +300x250 & 32 & 50x50 +300x250 & 32 & 100x75 diff --git a/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=5.output b/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=5.output new file mode 100644 index 0000000000..ec9538e8ac --- /dev/null +++ b/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=5.output @@ -0,0 +1,72 @@ +100x75 & 1 & 1x1 +100x75 & 1 & 10x10 +100x75 & 1 & 50x50 +100x75 & 1 & 100x75 +100x75 & 16 & 1x1 +100x75 & 16 & 10x10 +100x75 & 16 & 50x50 +100x75 & 16 & 100x75 +100x75 & 32 & 1x1 +100x75 & 32 & 10x10 +100x75 & 32 & 50x50 +100x75 & 32 & 100x75 +200x225 & 1 & 1x1 +200x225 & 1 & 10x10 +200x225 & 1 & 50x50 +200x225 & 1 & 100x75 +200x225 & 16 & 1x1 +200x225 & 16 & 10x10 +200x225 & 16 & 50x50 +200x225 & 16 & 100x75 +200x225 & 32 & 1x1 +200x225 & 32 & 10x10 +200x225 & 32 & 50x50 +200x225 & 32 & 100x75 +300x250 & 1 & 1x1 +300x250 & 1 & 10x10 +300x250 & 1 & 50x50 +300x250 & 1 & 100x75 +300x250 & 16 & 1x1 +300x250 & 16 & 10x10 +300x250 & 16 & 50x50 +300x250 & 16 & 100x75 +300x250 & 32 & 1x1 +300x250 & 32 & 10x10 +300x250 & 32 & 50x50 +300x250 & 32 & 100x75 +100x75 & 1 & 1x1 +100x75 & 1 & 10x10 +100x75 & 1 & 50x50 +100x75 & 1 & 100x75 +100x75 & 16 & 1x1 +100x75 & 16 & 10x10 +100x75 & 16 & 50x50 +100x75 & 16 & 100x75 +100x75 & 32 & 1x1 +100x75 & 32 & 10x10 +100x75 & 32 & 50x50 +100x75 & 32 & 100x75 +200x225 & 1 & 1x1 +200x225 & 1 & 10x10 +200x225 & 1 & 50x50 +200x225 & 1 & 100x75 +200x225 & 16 & 1x1 +200x225 & 16 & 10x10 +200x225 & 16 & 50x50 +200x225 & 16 & 100x75 +200x225 & 32 & 1x1 +200x225 & 32 & 10x10 +200x225 & 32 & 50x50 +200x225 & 32 & 100x75 +300x250 & 1 & 1x1 +300x250 & 1 & 10x10 +300x250 & 1 & 50x50 +300x250 & 1 & 100x75 +300x250 & 16 & 1x1 +300x250 & 16 & 10x10 +300x250 & 16 & 50x50 +300x250 & 16 & 100x75 +300x250 & 32 & 1x1 +300x250 & 32 & 10x10 +300x250 & 32 & 50x50 +300x250 & 32 & 100x75 -- 2.39.5