From 4b46c151c0c00e30f2842a2470c3ebcd586e79be Mon Sep 17 00:00:00 2001 From: Nicolas Barnafi Date: Wed, 7 Sep 2022 09:51:48 +0200 Subject: [PATCH] Fixed bddc test, added points as well --- include/deal.II/lac/petsc_precondition.h | 14 ++++---- source/lac/petsc_precondition.cc | 3 +- tests/petsc/bddc.cc | 41 ++++++++++++++++-------- 3 files changed, 37 insertions(+), 21 deletions(-) diff --git a/include/deal.II/lac/petsc_precondition.h b/include/deal.II/lac/petsc_precondition.h index 70416e60f0..6ac0b8a1f9 100644 --- a/include/deal.II/lac/petsc_precondition.h +++ b/include/deal.II/lac/petsc_precondition.h @@ -989,11 +989,11 @@ namespace PETScWrappers * Constructor. Note that BDDC offers a lot more options to set * 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, - std::vector> coords = {}); + AdditionalData(const bool use_vertices = true, + const bool use_edges = false, + const bool use_faces = false, + const bool symmetric = false, + const std::vector> coords = {}); /** * This flag sets the use of degrees of freedom in the vertices of the @@ -1024,7 +1024,7 @@ namespace PETScWrappers * Set the location of each DoF. This helps in improving the definition of * the vertices for unstructured meshes. */ - const std::vector> coords; + std::vector> coords; }; /** @@ -1079,6 +1079,8 @@ namespace PETScWrappers using PreconditionerBase DEAL_II_DEPRECATED = PreconditionBase; } // namespace PETScWrappers +template class PETScWrappers::PreconditionBDDC<2>; +template class PETScWrappers::PreconditionBDDC<3>; DEAL_II_NAMESPACE_CLOSE diff --git a/source/lac/petsc_precondition.cc b/source/lac/petsc_precondition.cc index 7950579cc6..467fa4ed43 100644 --- a/source/lac/petsc_precondition.cc +++ b/source/lac/petsc_precondition.cc @@ -912,8 +912,7 @@ namespace PETScWrappers set_option_value("-pc_bddc_symmetric", "true"); else set_option_value("-pc_bddc_symmetric", "false"); - if (additional_data.coords_cdim && additional_data.coords_n && - additional_data.coords_data) + if (additional_data.coords.size() > 0) { set_option_value("-pc_bddc_corner_selection", "true"); // Convert coords vector to PETSc data array diff --git a/tests/petsc/bddc.cc b/tests/petsc/bddc.cc index b5778ad23e..1dfc9fb76e 100644 --- a/tests/petsc/bddc.cc +++ b/tests/petsc/bddc.cc @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 2016 - 2021 by the deal.II authors +// Copyright (C) 2016 - 2022 by the deal.II authors // // This file is part of the deal.II library. // @@ -13,9 +13,10 @@ // // --------------------------------------------------------------------- -// Test that it is possible to instantiate a LinearOperator object for all -// different kinds of PETSc matrices and vectors -// TODO: A bit more tests... +// This test provides a minimal Laplace equation in 2D for testing +// the BDDC preconditioner. It cannot be simpler as it fails to setup +// in simpler matrices (i.e. diagonal), and the preset matrix factories +// create non-IS matrices, which are not supported by BDDC. #include @@ -44,6 +45,7 @@ #include #include +#include #include @@ -87,12 +89,12 @@ main(int argc, char *argv[]) locally_relevant_dofs); PETScWrappers::MPI::SparseMatrix system_matrix; - system_matrix.reinit_IS(locally_owned_dofs, - locally_active_dofs, - locally_owned_dofs, - locally_active_dofs, - dsp, - MPI_COMM_WORLD); + system_matrix.reinit(locally_owned_dofs, + locally_active_dofs, + locally_owned_dofs, + locally_active_dofs, + dsp, + MPI_COMM_WORLD); deallog << "MATIS:OK" << std::endl; PETScWrappers::MPI::Vector locally_relevant_solution(locally_owned_dofs, locally_relevant_dofs, @@ -154,9 +156,22 @@ main(int argc, char *argv[]) SolverControl solver_control(dof_handler.n_dofs(), 1e-12); - PETScWrappers::SolverCG solver(solver_control, MPI_COMM_WORLD); - PETScWrappers::PreconditionBDDC preconditioner(system_matrix); - + PETScWrappers::SolverCG solver(solver_control, MPI_COMM_WORLD); + PETScWrappers::PreconditionBDDC<2> preconditioner; + PETScWrappers::PreconditionBDDC<2>::AdditionalData data; + + // Now we setup the dof coordinates + 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) + { + coords[d2p.first] = d2p.second; + } + + preconditioner.initialize(system_matrix, data); check_solver_within_range(solver.solve(system_matrix, completely_distributed_solution, system_rhs, -- 2.39.5