From 8729b1939c42c0a92a8f5f8cc063abcfe220becf Mon Sep 17 00:00:00 2001
From: Nicolas Barnafi <nabw91@gmail.com>
Date: Thu, 8 Sep 2022 11:22:08 +0200
Subject: [PATCH] Added a check on the matrix type and a citation to the
 relevant reference.

---
 doc/doxygen/references.bib                    |  9 +++++++++
 doc/news/changes/minor/20221005NicolasBarnafi |  2 +-
 include/deal.II/lac/petsc_precondition.h      |  4 +---
 source/lac/petsc_precondition.cc              | 16 ++++++++++++++++
 4 files changed, 27 insertions(+), 4 deletions(-)

diff --git a/doc/doxygen/references.bib b/doc/doxygen/references.bib
index 136456e62f..69f1336238 100644
--- a/doc/doxygen/references.bib
+++ b/doc/doxygen/references.bib
@@ -1888,3 +1888,12 @@ url = {https://doi.org/10.1016/0045-7930(73)90027-3}
   doi       = {10.1007/978-3-031-07312-0_7},
   url       = {https://doi.org/10.1007/978-3-031-07312-0_7}
 }
+
+@misc{Barnafi2022,
+  doi = {10.48550/ARXIV.2208.06191},
+  url = {https://arxiv.org/abs/2208.06191},
+  author = {Barnafi, Nicolás A. and Pavarino, Luca F. and Scacchi, Simone},
+  title = {A comparative study of scalable multilevel preconditioners for cardiac mechanics},
+  publisher = {arXiv},
+  year = {2022},
+}
diff --git a/doc/news/changes/minor/20221005NicolasBarnafi b/doc/news/changes/minor/20221005NicolasBarnafi
index 6c930d10ed..e7e4cd151d 100644
--- a/doc/news/changes/minor/20221005NicolasBarnafi
+++ b/doc/news/changes/minor/20221005NicolasBarnafi
@@ -1,3 +1,3 @@
-New: The PetscWrappers::PreconditionBDDC preconditioner has been implemented. It requires a special matrix type (IS), which is initialized with the new function PetscWrappers::reinit_IS.
+New: The PetscWrappers::PreconditionBDDC preconditioner has been implemented. It requires a special matrix type (IS), which is initialized with the new function PetscWrappers::reinit.
 <br>
 (Nicolas Barnafi, 2022/05/10)
diff --git a/include/deal.II/lac/petsc_precondition.h b/include/deal.II/lac/petsc_precondition.h
index 6ac0b8a1f9..a4fa95a2ec 100644
--- a/include/deal.II/lac/petsc_precondition.h
+++ b/include/deal.II/lac/petsc_precondition.h
@@ -971,7 +971,7 @@ namespace PETScWrappers
    *   <li> A coarse solver: Continuity between each subdomain is imposed in a small number of DoFs, referred to as <em>primal DoFs</em>. This solver solves such problem.
    * </ul>
    *
-   * The size of the primal space is determined through the @p AdditionalData parameters.
+   * The size of the primal space is determined through the @p AdditionalData parameters. A thorough study of the performance of this solver in the context of cardiac mechanics, together with further details on this interface, is available in @cite Barnafi2022.
    *
    * @ingroup PETScWrappers
    */
@@ -1079,8 +1079,6 @@ 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 467fa4ed43..8c72ccd42e 100644
--- a/source/lac/petsc_precondition.cc
+++ b/source/lac/petsc_precondition.cc
@@ -894,6 +894,19 @@ namespace PETScWrappers
     PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCBDDC));
     AssertThrow(ierr == 0, ExcPETScError(ierr));
 
+    // The matrix must be of IS type. We check for this to avoid the PETSc error
+    // in order to suggest the correct matrix reinit method.
+    {
+      MatType current_type;
+      ierr = MatGetType(matrix, &current_type);
+      AssertThrow(ierr == 0, ExcPETScError(ierr));
+      AssertThrow(
+        strcmp(current_type, MATIS) == 0,
+        ExcMessage(
+          "Matrix must be of IS type. For this, the variant of reinit that includes the active dofs must be used."));
+    }
+
+
     std::stringstream ssStream;
 
     if (additional_data.use_vertices)
@@ -962,6 +975,9 @@ namespace PETScWrappers
 
 } // namespace PETScWrappers
 
+template class PETScWrappers::PreconditionBDDC<2>;
+template class PETScWrappers::PreconditionBDDC<3>;
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif // DEAL_II_WITH_PETSC
-- 
2.39.5