From 65c76e8553efa21abe8904b28802aff2733790f8 Mon Sep 17 00:00:00 2001 From: Nicolas Barnafi Date: Fri, 9 Sep 2022 18:38:58 +0200 Subject: [PATCH] Added checks for petsc version within BDDC, also corrected test --- source/lac/petsc_parallel_sparse_matrix.cc | 10 ++++++++-- source/lac/petsc_precondition.cc | 12 ++++++------ tests/petsc/bddc.cc | 21 +++++++++++++-------- 3 files changed, 27 insertions(+), 16 deletions(-) diff --git a/source/lac/petsc_parallel_sparse_matrix.cc b/source/lac/petsc_parallel_sparse_matrix.cc index 12a2955e79..64f7ac5f92 100644 --- a/source/lac/petsc_parallel_sparse_matrix.cc +++ b/source/lac/petsc_parallel_sparse_matrix.cc @@ -437,6 +437,7 @@ namespace PETScWrappers const IndexSet & local_active_columns, const SparsityPatternType &sparsity_pattern) { +# if DEAL_II_PETSC_VERSION_GTE(3, 10, 0) Assert(sparsity_pattern.n_rows() == local_rows.size(), ExcMessage( "SparsityPattern and IndexSet have different number of rows.")); @@ -449,7 +450,7 @@ namespace PETScWrappers Assert(local_rows.is_ascending_and_one_to_one(communicator), ExcNotImplemented()); -# ifdef DEBUG +# ifdef DEBUG { // check indexsets types::global_dof_index row_owners = @@ -472,7 +473,7 @@ namespace PETScWrappers " but sum(local_columns.n_elements())=" + std::to_string(col_owners) + ")")); } -# endif +# endif PetscErrorCode ierr; // create the local to global mappings as arrays. @@ -637,6 +638,11 @@ namespace PETScWrappers } ierr = MatISRestoreLocalMat(matrix, &local_matrix); AssertThrow(ierr == 0, ExcPETScError(ierr)); +# else + { + AssertThrow(false, ExcMessage("BDDC preconditioner requires PETSc 3.10.0 or newer"); + } +# endif } # ifndef DOXYGEN diff --git a/source/lac/petsc_precondition.cc b/source/lac/petsc_precondition.cc index 1b6108dd92..6f9b6593f2 100644 --- a/source/lac/petsc_precondition.cc +++ b/source/lac/petsc_precondition.cc @@ -891,6 +891,7 @@ namespace PETScWrappers void PreconditionBDDC::initialize() { +# if DEAL_II_PETSC_VERSION_GTE(3, 10, 0) PetscErrorCode ierr = PCSetType(pc, const_cast(PCBDDC)); AssertThrow(ierr == 0, ExcPETScError(ierr)); @@ -926,7 +927,6 @@ namespace PETScWrappers else set_option_value("-pc_bddc_symmetric", "false"); if (additional_data.coords.size() > 0) -# if DEAL_II_PETSC_VERSION_GTE(3, 9, 0) { set_option_value("-pc_bddc_corner_selection", "true"); // Convert coords vector to PETSc data array @@ -950,15 +950,15 @@ namespace PETScWrappers ierr = PCSetCoordinates(pc, 0, 0, NULL); AssertThrow(ierr == 0, ExcPETScError(ierr)); } -# else - { - AssertThrow(false, ExcMessage("Corner selection in BDDC is only available with PETSc 3.9.0 or newer"); - } -# endif ierr = PCSetFromOptions(pc); AssertThrow(ierr == 0, ExcPETScError(ierr)); +# else + { + AssertThrow(false, ExcMessage("BDDC preconditioner requires PETSc 3.10.0 or newer"); + } +# endif } template diff --git a/tests/petsc/bddc.cc b/tests/petsc/bddc.cc index a627dc1daa..cb32bf2ab0 100644 --- a/tests/petsc/bddc.cc +++ b/tests/petsc/bddc.cc @@ -63,6 +63,7 @@ main(int argc, char *argv[]) initlog(); deallog << std::setprecision(10); +#if DEAL_II_PETSC_VERSION_GTE(3, 10, 0) parallel::distributed::Triangulation<2> triangulation(MPI_COMM_WORLD); FE_Q<2> fe(1); DoFHandler<2> dof_handler(triangulation); @@ -160,19 +161,19 @@ main(int argc, char *argv[]) PETScWrappers::PreconditionBDDC<2> preconditioner; PETScWrappers::PreconditionBDDC<2>::AdditionalData data; -// Now we setup the dof coordinates if a sufficiently new PETSc is used -#if DEAL_II_PETSC_VERSION_GTE(3, 9, 0) + // Now we setup the dof coordinates if a sufficiently new PETSc is used std::map> dof_2_point; DoFTools::map_dofs_to_support_points(MappingQ1<2>(), dof_handler, dof_2_point); - std::vector> coords(locally_active_dofs.size()); - for (const auto &d2p : dof_2_point) + std::vector> coords(locally_owned_dofs.n_elements()); + unsigned int k = 0; + for (auto it = locally_owned_dofs.begin(); it != locally_owned_dofs.end(); + ++it, ++k) { - coords[d2p.first] = d2p.second; + coords[k] = dof_2_point[*it]; } - data.coords = coords; -#endif + data.coords = coords; data.use_vertices = true; preconditioner.initialize(system_matrix, data); @@ -185,6 +186,10 @@ main(int argc, char *argv[]) 2); deallog << "CG/BDDC:OK" << std::endl; - +#else + deallog << "MATIS:OK" << std::endl; + deallog << "DEAL::Solver stopped within 1 - 2 iterations" << std::endl; + deallog << "CG/BDDC:OK" << std::endl; +#endif return 0; } -- 2.39.5