From 14b3e66fa59045443b3a4011b5990ea379f2bae0 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Tue, 9 Jan 2024 15:18:58 -0700 Subject: [PATCH] Minor doc updates to the sparse direct solvers. --- include/deal.II/lac/petsc_solver.h | 26 ++++++++++++++++++++++---- include/deal.II/lac/sparse_direct.h | 5 +++++ source/lac/petsc_solver.cc | 4 ++-- 3 files changed, 29 insertions(+), 6 deletions(-) diff --git a/include/deal.II/lac/petsc_solver.h b/include/deal.II/lac/petsc_solver.h index 51ef469f69..ac95810089 100644 --- a/include/deal.II/lac/petsc_solver.h +++ b/include/deal.II/lac/petsc_solver.h @@ -934,12 +934,30 @@ namespace PETScWrappers solve(const MatrixBase &A, VectorBase &x, const VectorBase &b); /** - * The method allows to take advantage if the system matrix is symmetric - * by using LDL^T decomposition instead of more expensive LU. The argument - * indicates whether the matrix is symmetric or not. + * If called with `true` as argument, tell the direct solver + * to assume that the system matrix is symmetric. + * It does so by computing the LDL^T decomposition (in effect, a + * Cholesky decomposition) instead of more expensive LU decomposition. + * The argument indicates whether the matrix can be assumed to be + * symmetric or not. + * + * Note that most finite element matrices are "structurally symmetric", + * i.e., the *sparsity pattern* is symmetric, even though the matrix + * is not. An example of a matrix that is structurally symmetric + * but not symmetric is the matrix you obtain by discretizing the + * advection equation $\nabla \cdot (\vec\beta u) = f$ (see, for + * example step-12). Because the operator here is not symmetric, the + * matrix is not symmetric either; however, if matrix entry $A_{ij}$ + * is nonzero, then matrix entry $A_{ji}$ is generally not zero either + * (and in any case, DoFTools::make_sparsity_pattern() will create a + * symmetric sparsity pattern). That said, the current function + * is not meant to indicate whether the *sparsity* pattern is symmetric, + * but whether the *matrix* itself is symmetric, and this typically + * requires that the differential operator you are considering is + * symmetric. */ void - set_symmetric_mode(const bool flag); + set_symmetric_mode(const bool matrix_is_symmetric); protected: /** diff --git a/include/deal.II/lac/sparse_direct.h b/include/deal.II/lac/sparse_direct.h index 8cc8e200c7..ee3377723a 100644 --- a/include/deal.II/lac/sparse_direct.h +++ b/include/deal.II/lac/sparse_direct.h @@ -82,6 +82,11 @@ namespace types * SparseMatrix, SparseMatrixEZ, SparseMatrixEZ, * BlockSparseMatrix, and BlockSparseMatrix. * + * This class is not instantiated for the matrix types in namespace + * PETScWrappers or TrilinosWrappers. However, PETScWrappers::SparseDirectMUMPS + * can be used for PETSc matrices, and in fact even works for parallel + * computations. + * * @ingroup Solvers Preconditioners */ class SparseDirectUMFPACK : public Subscriptor diff --git a/source/lac/petsc_solver.cc b/source/lac/petsc_solver.cc index 3944a275cf..857dc9b951 100644 --- a/source/lac/petsc_solver.cc +++ b/source/lac/petsc_solver.cc @@ -812,9 +812,9 @@ namespace PETScWrappers void - SparseDirectMUMPS::set_symmetric_mode(const bool flag) + SparseDirectMUMPS::set_symmetric_mode(const bool matrix_is_symmetric) { - symmetric_mode = flag; + symmetric_mode = matrix_is_symmetric; } } // namespace PETScWrappers -- 2.39.5