]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Minor doc updates to the sparse direct solvers. 16443/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Tue, 9 Jan 2024 22:18:58 +0000 (15:18 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Tue, 9 Jan 2024 22:18:58 +0000 (15:18 -0700)
include/deal.II/lac/petsc_solver.h
include/deal.II/lac/sparse_direct.h
source/lac/petsc_solver.cc

index 51ef469f69f98e4b34b7b1e187276957ab752df3..ac95810089ea9d195e1c72ca4ef330fba0825c76 100644 (file)
@@ -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:
     /**
index 8cc8e200c7b40879e988731a1a3963cc882841d9..ee3377723a2431c434d23e64a71c23726ed04ef9 100644 (file)
@@ -82,6 +82,11 @@ namespace types
  * SparseMatrix<float>, SparseMatrixEZ<float>, SparseMatrixEZ<double>,
  * BlockSparseMatrix<double>, and BlockSparseMatrix<float>.
  *
+ * 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
index 3944a275cfa037b4bebda55050b0bcd68c388c0c..857dc9b951ae49a1d65ea706c3c56f9847d21435 100644 (file)
@@ -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

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.