From 5779d3263686039ce15a3758b51aa0b05116dbc4 Mon Sep 17 00:00:00 2001 From: Alexander Grayver Date: Fri, 3 Oct 2014 19:57:50 +0200 Subject: [PATCH] Enable PETScWrappers::MPI::SparseMatrix interface to account for off-diagonal elements during preallocation. --- .../lac/petsc_parallel_sparse_matrix.h | 38 ++++++++++------ source/lac/petsc_parallel_sparse_matrix.cc | 43 +++++++++++++------ 2 files changed, 54 insertions(+), 27 deletions(-) diff --git a/include/deal.II/lac/petsc_parallel_sparse_matrix.h b/include/deal.II/lac/petsc_parallel_sparse_matrix.h index a28dc035ff..1bf0bc01cd 100644 --- a/include/deal.II/lac/petsc_parallel_sparse_matrix.h +++ b/include/deal.II/lac/petsc_parallel_sparse_matrix.h @@ -166,11 +166,15 @@ namespace PETScWrappers * Create a sparse matrix of * dimensions @p m times @p n, with * an initial guess of @p - * n_nonzero_per_row nonzero elements - * per row. PETSc is able to cope - * with the situation that more than - * this number of elements are later - * allocated for a row, but this + * n_nonzero_per_row and @p + * n_offdiag_nonzero_per_row nonzero + * elements per row (see documentation + * of the MatCreateAIJ PETSc function + * for more information about these + * parameters). PETSc is able to + * cope with the situation that more + * than this number of elements are + * later allocated for a row, but this * involves copying data, and is thus * expensive. * @@ -198,14 +202,17 @@ namespace PETScWrappers const size_type local_rows, const size_type local_columns, const size_type n_nonzero_per_row, - const bool is_symmetric = false); + const bool is_symmetric = false, + const size_type n_offdiag_nonzero_per_row = 0); /** * Initialize a rectangular matrix * with @p m rows and @p n columns. * The maximal number of nonzero - * entries for each row separately is - * given by the @p row_lengths array. + * entries for diagonal and off- + * diagonal blocks of each row is + * given by the @p row_lengths and + * @p offdiag_row_lengths arrays. * * For the meaning of the @p * local_row and @p local_columns @@ -239,7 +246,8 @@ namespace PETScWrappers const size_type local_rows, const size_type local_columns, const std::vector &row_lengths, - const bool is_symmetric = false); + const bool is_symmetric = false, + const std::vector &offdiag_row_lengths = std::vector()); /** * Initialize using the given @@ -328,7 +336,8 @@ namespace PETScWrappers const size_type local_rows, const size_type local_columns, const size_type n_nonzero_per_row, - const bool is_symmetric = false); + const bool is_symmetric = false, + const size_type n_offdiag_nonzero_per_row = 0); /** * Throw away the present matrix and @@ -344,7 +353,8 @@ namespace PETScWrappers const size_type local_rows, const size_type local_columns, const std::vector &row_lengths, - const bool is_symmetric = false); + const bool is_symmetric = false, + const std::vector &offdiag_row_lengths = std::vector()); /** * Initialize using the given @@ -482,7 +492,8 @@ namespace PETScWrappers const size_type local_rows, const size_type local_columns, const size_type n_nonzero_per_row, - const bool is_symmetric = false); + const bool is_symmetric = false, + const size_type n_offdiag_nonzero_per_row = 0); /** * Same as previous function. @@ -492,7 +503,8 @@ namespace PETScWrappers const size_type local_rows, const size_type local_columns, const std::vector &row_lengths, - const bool is_symmetric = false); + const bool is_symmetric = false, + const std::vector &offdiag_row_lengths = std::vector()); /** * Same as previous functions. diff --git a/source/lac/petsc_parallel_sparse_matrix.cc b/source/lac/petsc_parallel_sparse_matrix.cc index 7a842188a0..3e59ffe75e 100644 --- a/source/lac/petsc_parallel_sparse_matrix.cc +++ b/source/lac/petsc_parallel_sparse_matrix.cc @@ -60,12 +60,14 @@ namespace PETScWrappers const size_type local_rows, const size_type local_columns, const size_type n_nonzero_per_row, - const bool is_symmetric) + const bool is_symmetric, + const size_type n_offdiag_nonzero_per_row) : communicator (communicator) { do_reinit (m, n, local_rows, local_columns, - n_nonzero_per_row, is_symmetric); + n_nonzero_per_row, is_symmetric, + n_offdiag_nonzero_per_row); } @@ -76,12 +78,13 @@ namespace PETScWrappers const size_type local_rows, const size_type local_columns, const std::vector &row_lengths, - const bool is_symmetric) + const bool is_symmetric, + const std::vector &offdiag_row_lengths) : communicator (communicator) { do_reinit (m, n, local_rows, local_columns, - row_lengths, is_symmetric); + row_lengths, is_symmetric, offdiag_row_lengths); } @@ -130,7 +133,8 @@ namespace PETScWrappers const size_type local_rows, const size_type local_columns, const size_type n_nonzero_per_row, - const bool is_symmetric) + const bool is_symmetric, + const size_type n_offdiag_nonzero_per_row) { this->communicator = communicator; @@ -144,7 +148,8 @@ namespace PETScWrappers AssertThrow (ierr == 0, ExcPETScError(ierr)); do_reinit (m, n, local_rows, local_columns, - n_nonzero_per_row, is_symmetric); + n_nonzero_per_row, is_symmetric, + n_offdiag_nonzero_per_row); } @@ -156,7 +161,8 @@ namespace PETScWrappers const size_type local_rows, const size_type local_columns, const std::vector &row_lengths, - const bool is_symmetric) + const bool is_symmetric, + const std::vector &offdiag_row_lengths) { this->communicator = communicator; @@ -169,7 +175,8 @@ namespace PETScWrappers #endif AssertThrow (ierr == 0, ExcPETScError(ierr)); - do_reinit (m, n, local_rows, local_columns, row_lengths, is_symmetric); + do_reinit (m, n, local_rows, local_columns, + row_lengths, is_symmetric, offdiag_row_lengths); } @@ -228,7 +235,8 @@ namespace PETScWrappers const size_type local_rows, const size_type local_columns, const size_type n_nonzero_per_row, - const bool is_symmetric) + const bool is_symmetric, + const size_type n_offdiag_nonzero_per_row) { Assert (local_rows <= m, ExcLocalRowsTooLarge (local_rows, m)); @@ -242,14 +250,16 @@ namespace PETScWrappers = MatCreateMPIAIJ (communicator, local_rows, local_columns, m, n, - n_nonzero_per_row, 0, 0, 0, + n_nonzero_per_row, 0, + n_offdiag_nonzero_per_row, 0, &matrix); #else ierr = MatCreateAIJ (communicator, local_rows, local_columns, m, n, - n_nonzero_per_row, 0, 0, 0, + n_nonzero_per_row, 0, + n_offdiag_nonzero_per_row, 0, &matrix); AssertThrow (ierr == 0, ExcPETScError(ierr)); @@ -280,7 +290,8 @@ namespace PETScWrappers const size_type local_rows, const size_type local_columns, const std::vector &row_lengths, - const bool is_symmetric) + const bool is_symmetric, + const std::vector &offdiag_row_lengths) { Assert (local_rows <= m, ExcLocalRowsTooLarge (local_rows, m)); @@ -307,6 +318,8 @@ namespace PETScWrappers // tricks with conversions of pointers const std::vector int_row_lengths (row_lengths.begin(), row_lengths.end()); + const std::vector int_offdiag_row_lengths (offdiag_row_lengths.begin(), + offdiag_row_lengths.end()); //TODO: There must be a significantly better way to provide information about the off-diagonal blocks of the matrix. this way, petsc keeps allocating tiny chunks of memory, and gets completely hung up over this @@ -317,14 +330,16 @@ namespace PETScWrappers = MatCreateMPIAIJ (communicator, local_rows, local_columns, m, n, - 0, &int_row_lengths[0], 0, 0, + 0, &int_row_lengths[0], + 0, offdiag_row_lengths.size() ? &int_offdiag_row_lengths[0] : 0, &matrix); #else ierr = MatCreateAIJ (communicator, local_rows, local_columns, m, n, - 0, &int_row_lengths[0], 0, 0, + 0, &int_row_lengths[0], + 0, offdiag_row_lengths.size() ? &int_offdiag_row_lengths[0] : 0, &matrix); AssertThrow (ierr == 0, ExcPETScError(ierr)); -- 2.39.5