]> https://gitweb.dealii.org/ - dealii.git/commitdiff
address comments, add documentation
authorTimo Heister <timo.heister@gmail.com>
Mon, 23 Mar 2020 12:39:22 +0000 (08:39 -0400)
committerTimo Heister <timo.heister@gmail.com>
Mon, 23 Mar 2020 17:13:24 +0000 (13:13 -0400)
include/deal.II/lac/precondition.h

index 7f07e70f9fbbaac26a525472397dab40008b1c43..f2261ddc52c06abb07f92b153a863fab31f00bde 100644 (file)
@@ -914,7 +914,7 @@ private:
  * which are computed during the first invocation of vmult(). The algorithm
  * invokes a conjugate gradient solver (i.e., Lanczos iteration) so symmetry
  * and positive definiteness of the (preconditioned) matrix system are
- * requirements. The eigenvalue algorithm can be controlled by
+ * required. The eigenvalue algorithm can be controlled by
  * PreconditionChebyshev::AdditionalData::eig_cg_n_iterations specifying how
  * many iterations should be performed. The iterations are started from an
  * initial vector that depends on the vector type. For the classes
@@ -929,16 +929,20 @@ private:
  * length of the vector, except for the very first entry that is zero,
  * triggering high-frequency content again.
  *
- * The computation of eigenvalues happens the first time one of the
- * vmult(), Tvmult(), step() or Tstep() functions is called. This is because
- * temporary vectors of the same layout as the source and destination vectors
- * are necessary for these computations and this information gets only
- * available through vmult().
+ * The computation of eigenvalues happens the first time one of the vmult(),
+ * Tvmult(), step() or Tstep() functions is called or when
+ * estimate_eigenvalues() is called directly. In the latter case, it is
+ * necessary to provide a temporary vector of the same layout as the source
+ * and destination vectors used during application of the preconditioner.
  *
- * Due to the cost of the eigenvalue estimate in the first vmult(), this class
- * is most appropriate if it is applied repeatedly, e.g. in a smoother for a
- * geometric multigrid solver, that can in turn be used to solve several
- * linear systems.
+ * The estimates for minimum and maximum eigenvalue are taken from SolverCG
+ * (even if the solver did not converge in the requested number of
+ * iterations). Finally, the maximum eigenvalue is multiplied by a safety
+ * factor of 1.2.
+ *
+ * Due to the cost of the eigenvalue estimate, this class is most appropriate
+ * if it is applied repeatedly, e.g. in a smoother for a geometric multigrid
+ * solver, that can in turn be used to solve several linear systems.
  *
  * <h4>Bypassing the eigenvalue computation</h4>
  *
@@ -1161,25 +1165,26 @@ public:
   n() const;
 
   /**
-   * Struct that contains information about the eigenvalue estimation performed
-   * by this class.
+   * A struct that contains information about the eigenvalue estimation
+   * performed by the PreconditionChebychev class.
    */
   struct EigenvalueInformation
   {
     /**
      * Estimate for the smallest eigenvalue.
      */
-    double min_eigenvalue;
+    double min_eigenvalue_estimate;
     /**
      * Estimate for the largest eigenvalue.
      */
-    double max_eigenvalue;
+    double max_eigenvalue_estimate;
     /**
      * Number of CG iterations performed or 0.
      */
     unsigned int cg_iterations;
     /**
-     * The degree (either as set or estimated).
+     * The degree of the Chebyshev polynomial (either as set using
+     * AdditionalData::degree or estimated as described there).
      */
     unsigned int degree;
   };
@@ -2325,28 +2330,28 @@ PreconditionChebyshev<MatrixType, VectorType, PreconditionerType>::
 
       // read the eigenvalues from the attached eigenvalue tracker
       if (eigenvalue_tracker.values.empty())
-        info.min_eigenvalue = info.max_eigenvalue = 1.;
+        info.min_eigenvalue_estimate = info.max_eigenvalue_estimate = 1.;
       else
         {
-          info.min_eigenvalue = eigenvalue_tracker.values.front();
+          info.min_eigenvalue_estimate = eigenvalue_tracker.values.front();
 
           // include a safety factor since the CG method will in general not
           // be converged
-          info.max_eigenvalue = 1.2 * eigenvalue_tracker.values.back();
+          info.max_eigenvalue_estimate = 1.2 * eigenvalue_tracker.values.back();
         }
 
       info.cg_iterations = control.last_step();
     }
   else
     {
-      info.max_eigenvalue = data.max_eigenvalue;
-      info.min_eigenvalue = data.max_eigenvalue / data.smoothing_range;
+      info.max_eigenvalue_estimate = data.max_eigenvalue;
+      info.min_eigenvalue_estimate = data.max_eigenvalue / data.smoothing_range;
     }
 
-  const double alpha =
-    (data.smoothing_range > 1. ?
-       info.max_eigenvalue / data.smoothing_range :
-       std::min(0.9 * info.max_eigenvalue, info.min_eigenvalue));
+  const double alpha = (data.smoothing_range > 1. ?
+                          info.max_eigenvalue_estimate / data.smoothing_range :
+                          std::min(0.9 * info.max_eigenvalue_estimate,
+                                   info.min_eigenvalue_estimate));
 
   // in case the user set the degree to invalid unsigned int, we have to
   // determine the number of necessary iterations from the Chebyshev error
@@ -2355,7 +2360,7 @@ PreconditionChebyshev<MatrixType, VectorType, PreconditionerType>::
   // R. S. Varga, Matrix iterative analysis, 2nd ed., Springer, 2009
   if (data.degree == numbers::invalid_unsigned_int)
     {
-      const double actual_range = info.max_eigenvalue / alpha;
+      const double actual_range = info.max_eigenvalue_estimate / alpha;
       const double sigma        = (1. - std::sqrt(1. / actual_range)) /
                            (1. + std::sqrt(1. / actual_range));
       const double eps = data.smoothing_range;
@@ -2372,10 +2377,10 @@ PreconditionChebyshev<MatrixType, VectorType, PreconditionerType>::
 
   const_cast<
     PreconditionChebyshev<MatrixType, VectorType, PreconditionerType> *>(this)
-    ->delta = (info.max_eigenvalue - alpha) * 0.5;
+    ->delta = (info.max_eigenvalue_estimate - alpha) * 0.5;
   const_cast<
     PreconditionChebyshev<MatrixType, VectorType, PreconditionerType> *>(this)
-    ->theta = (info.max_eigenvalue + alpha) * 0.5;
+    ->theta = (info.max_eigenvalue_estimate + alpha) * 0.5;
 
   // We do not need the second temporary vector in case we have a
   // DiagonalMatrix as preconditioner and use deal.II's own vectors

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.