From 271310d9ebf57db3a3541a52feaf6cb4aad378fa Mon Sep 17 00:00:00 2001 From: Nicolas Barnafi Date: Mon, 4 Jul 2022 12:04:54 +0200 Subject: [PATCH] Added missing documentation, plus some missing ierr's and indentend --- include/deal.II/lac/petsc_precondition.h | 50 +++++++++---- source/lac/petsc_parallel_sparse_matrix.cc | 82 +++++++++++++--------- source/lac/petsc_precondition.cc | 4 +- 3 files changed, 87 insertions(+), 49 deletions(-) diff --git a/include/deal.II/lac/petsc_precondition.h b/include/deal.II/lac/petsc_precondition.h index 7a27046a8e..20a2c3a182 100644 --- a/include/deal.II/lac/petsc_precondition.h +++ b/include/deal.II/lac/petsc_precondition.h @@ -959,7 +959,18 @@ namespace PETScWrappers /** * A class that implements the interface to use the BDDC preconditioner from - * PETSc. + * PETSc (PCBDDC), + * which is a two-level, substructuring, non-overlapping domain decomposition + * preconditioner. Details of the implementation can be found in "S Zampini, + * SISC (2016)". It mainly consists of two elements: + * + * + * + * The size of the primal space is determined through the @p AdditionalData parameters. * * @ingroup PETScWrappers */ @@ -977,28 +988,30 @@ namespace PETScWrappers * than what is exposed here. */ AdditionalData(const bool use_vertices = true, - const bool use_edges_ = false, - const bool use_faces_ = false, - const bool symmetric_ = false, - const unsigned int coords_cdim_ = 0, - const types::global_dof_index coords_n_ = 0, - const PetscReal * coords_data_ = nullptr); + const bool use_edges_ = false, + const bool use_faces_ = false, + const bool symmetric_ = false, + const unsigned int coords_cdim_ = 0, + const types::global_dof_index coords_n_ = 0, + const PetscReal * coords_data_ = nullptr); /** * This flag sets the use of degrees of freedom in the vertices of the - * subdomain as primal variables. + * subdomains as primal variables for the creation of the coarse space. */ bool use_vertices; /** * This flag sets the use of degrees of freedom in the edges of the - * subdomain as primal variables. + * subdomain as primal variables for the creation of the coarse space. + * Continuity is actually imposed at the edge average. */ bool use_edges; /** * This flag sets the use of degrees of freedom in the faces of the - * subdomain as primal variables. + * subdomain as primal variables for the creation of the coarse space. + * Continuity is actually imposed at the face average. */ bool use_faces; @@ -1008,12 +1021,21 @@ namespace PETScWrappers bool symmetric; /** - * Support for H1 problems + * Set the number of coordinates that each DoF has, i.e. 2 in 2D and 3 in + * 3D. + */ + unsigned int coords_cdim; + + /** + * Set the number of degrees of freedom that the problem has. */ - unsigned int coords_cdim; types::global_dof_index coords_n; - const PetscReal - *coords_data; /* not referenced, consumed at initialize time */ + + /** + * Set the location of each DoF. This helps in improving the definition of + * the vertices for unstructured meshes. + */ + const PetscReal *coords_data; }; /** diff --git a/source/lac/petsc_parallel_sparse_matrix.cc b/source/lac/petsc_parallel_sparse_matrix.cc index fb7bc57667..fcc2576b6c 100644 --- a/source/lac/petsc_parallel_sparse_matrix.cc +++ b/source/lac/petsc_parallel_sparse_matrix.cc @@ -473,6 +473,7 @@ namespace PETScWrappers std::to_string(col_owners) + ")")); } # endif + PetscErrorCode ierr; // create the local to global mappings as arrays. IndexSet::size_type n_l2g_row = local_active_rows.n_elements(); @@ -493,39 +494,51 @@ namespace PETScWrappers IS is_glob_row, is_glob_col; // Create row index set ISLocalToGlobalMapping l2gmap_row; - ISCreateGeneral(communicator, - n_l2g_row, - idx_glob_row.data(), - PETSC_COPY_VALUES, - &is_glob_row); - ISLocalToGlobalMappingCreateIS(is_glob_row, &l2gmap_row); - ISDestroy(&is_glob_row); - ISLocalToGlobalMappingViewFromOptions(l2gmap_row, NULL, "-view_map"); + ierr = ISCreateGeneral(communicator, + n_l2g_row, + idx_glob_row.data(), + PETSC_COPY_VALUES, + &is_glob_row); + AssertThrow(ierr == 0, ExcPETScError(ierr)); + ierr = ISLocalToGlobalMappingCreateIS(is_glob_row, &l2gmap_row); + AssertThrow(ierr == 0, ExcPETScError(ierr)); + ierr = ISDestroy(&is_glob_row); + AssertThrow(ierr == 0, ExcPETScError(ierr)); + ierr = + ISLocalToGlobalMappingViewFromOptions(l2gmap_row, NULL, "-view_map"); + AssertThrow(ierr == 0, ExcPETScError(ierr)); // Create column index set ISLocalToGlobalMapping l2gmap_col; - ISCreateGeneral(communicator, - n_l2g_col, - idx_glob_col.data(), - PETSC_COPY_VALUES, - &is_glob_col); - ISLocalToGlobalMappingCreateIS(is_glob_col, &l2gmap_col); - ISDestroy(&is_glob_col); - ISLocalToGlobalMappingViewFromOptions(l2gmap_col, NULL, "-view_map"); + ierr = ISCreateGeneral(communicator, + n_l2g_col, + idx_glob_col.data(), + PETSC_COPY_VALUES, + &is_glob_col); + AssertThrow(ierr == 0, ExcPETScError(ierr)); + ierr = ISLocalToGlobalMappingCreateIS(is_glob_col, &l2gmap_col); + AssertThrow(ierr == 0, ExcPETScError(ierr)); + ierr = ISDestroy(&is_glob_col); + AssertThrow(ierr == 0, ExcPETScError(ierr)); + ierr = + ISLocalToGlobalMappingViewFromOptions(l2gmap_col, NULL, "-view_map"); + AssertThrow(ierr == 0, ExcPETScError(ierr)); // create the matrix with the IS constructor. - PetscErrorCode ierr = MatCreateIS(communicator, - 1, - local_rows.n_elements(), - local_columns.n_elements(), - sparsity_pattern.n_rows(), - sparsity_pattern.n_cols(), - l2gmap_row, - l2gmap_col, - &matrix); - AssertThrow(ierr == 0, ExcPETScError(ierr)); - ISLocalToGlobalMappingDestroy(&l2gmap_row); - ISLocalToGlobalMappingDestroy(&l2gmap_col); + ierr = MatCreateIS(communicator, + 1, + local_rows.n_elements(), + local_columns.n_elements(), + sparsity_pattern.n_rows(), + sparsity_pattern.n_cols(), + l2gmap_row, + l2gmap_col, + &matrix); + AssertThrow(ierr == 0, ExcPETScError(ierr)); + ierr = ISLocalToGlobalMappingDestroy(&l2gmap_row); + AssertThrow(ierr == 0, ExcPETScError(ierr)); + ierr = ISLocalToGlobalMappingDestroy(&l2gmap_col); + AssertThrow(ierr == 0, ExcPETScError(ierr)); // next preset the exact given matrix // entries with zeros. This doesn't avoid any @@ -542,10 +555,12 @@ namespace PETScWrappers // In the MATIS case, we use the local matrix instead Mat local_matrix; - MatISGetLocalMat(matrix, &local_matrix); - MatSetType(local_matrix, - MATSEQAIJ); // SEQ as it is local! TODO: Allow for OpenMP - // parallelization in local node. + ierr = MatISGetLocalMat(matrix, &local_matrix); + AssertThrow(ierr == 0, ExcPETScError(ierr)); + ierr = MatSetType(local_matrix, + MATSEQAIJ); // SEQ as it is local! TODO: Allow for + // OpenMP parallelization in local node. + AssertThrow(ierr == 0, ExcPETScError(ierr)); if (local_rows.n_elements() > 0) { // MatSEQAIJSetPreallocationCSR @@ -621,7 +636,8 @@ namespace PETScWrappers close_matrix(local_matrix); set_keep_zero_rows(local_matrix); } - MatISRestoreLocalMat(matrix, &local_matrix); + ierr = MatISRestoreLocalMat(matrix, &local_matrix); + AssertThrow(ierr == 0, ExcPETScError(ierr)); } # ifndef DOXYGEN diff --git a/source/lac/petsc_precondition.cc b/source/lac/petsc_precondition.cc index aa2cd5caf9..f78ece4248 100644 --- a/source/lac/petsc_precondition.cc +++ b/source/lac/petsc_precondition.cc @@ -892,7 +892,6 @@ namespace PETScWrappers PetscErrorCode ierr = PCSetType(pc, const_cast(PCBDDC)); AssertThrow(ierr == 0, ExcPETScError(ierr)); - // TODO: Add BDDC options, in particular the additional primal nodes std::stringstream ssStream; if (additional_data.use_vertices) @@ -911,7 +910,8 @@ namespace PETScWrappers set_option_value("-pc_bddc_symmetric", "true"); else set_option_value("-pc_bddc_symmetric", "false"); - if (additional_data.coords_cdim) + if (additional_data.coords_cdim && additional_data.coords_n && + additional_data.coords_data) { set_option_value("-pc_bddc_corner_selection", "true"); ierr = PCSetCoordinates(pc, -- 2.39.5