]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Augment documentation
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Fri, 9 Jun 2017 06:19:01 +0000 (08:19 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Fri, 9 Jun 2017 16:19:38 +0000 (18:19 +0200)
include/deal.II/multigrid/mg_base.h
include/deal.II/multigrid/mg_smoother.h
include/deal.II/multigrid/multigrid.h

index d825ed8b29b4e6b5a4318a2560793d05055d4bf9..d1973a4a25abe1c01559a0040aaec80953b24a11 100644 (file)
@@ -213,6 +213,21 @@ public:
  * Base class for multigrid smoothers. Does nothing but defining the interface
  * used by multigrid methods.
  *
+ * The smoother interface provides two methods, a smooth() method and an
+ * apply() method. The multigrid preconditioner interfaces distinguish between
+ * the two for efficiency reasons: Upon entry to the preconditioner operations,
+ * the vector @p u needs to be set to zero and smoothing starts by a simple
+ * application of the smoother on the @p rhs vector. This method is provided by
+ * the apply() method of this class. It is the same as first setting @p u to
+ * zero and then calling smooth(), but for many classes the separate apply()
+ * interface is more efficient because it can skip one matrix-vector product.
+ *
+ * In the multigrid preconditioner interfaces, the apply() method is used for
+ * the pre-smoothing operation because the previous content in the solution
+ * vector needs to be overwritten for a new incoming residual. On the other
+ * hand, all subsequent operations need to smooth the content already present
+ * in the vector @p u given the right hand side, which is done by smooth().
+ *
  * @author Guido Kanschat, 2002
  */
 template <typename VectorType>
@@ -230,21 +245,30 @@ public:
   virtual void clear() = 0;
 
   /**
-   * Smoothing function. This is the function used in multigrid methods.
+   * Smoothing function that smooths the content in @p u given the right hand
+   * side vector @p rhs. This is the function used in multigrid methods.
    */
   virtual void smooth (const unsigned int level,
                        VectorType         &u,
                        const VectorType   &rhs) const = 0;
 
   /**
-   * As opposed to the smoothing function, this function applies the action of
+   * As opposed to the smooth() function, this function applies the action of
    * the smoothing, overwriting the previous content in the vector u. This
    * function must be equivalent to the following code
    * @code
    * u = 0;
    * smooth(level, u, rhs);
    * @endcode
-   * but can usually be implemented more efficiently than the former.
+   * but can usually be implemented more efficiently than the former. If a
+   * particular smoother does not override the apply() method, the default
+   * implementation as described here is used.
+   *
+   * In the multigrid preconditioner interfaces, the apply() method is used for
+   * the pre-smoothing operation because the previous content in the solution
+   * vector needs to be overwritten for a new incoming residual. On the other
+   * hand, all subsequent operations need to smooth the content already present
+   * in the vector @p u given the right hand side, which is done by smooth().
    */
   virtual void apply (const unsigned int level,
                       VectorType         &u,
index 11f4236927bfb20da8e5a07cc5653e1fb105f6e8..0c00827513ee17457ada4f04d02dd7091dc9ca2d 100644 (file)
@@ -227,7 +227,18 @@ namespace mg
 
     /**
      * The apply variant of smoothing, setting the vector u to zero before
-     * calling the smooth function.
+     * calling the smooth function. This function is equivalent to the
+     * following code
+     * @code
+     * u = 0;
+     * smooth(level, u, rhs);
+     * @endcode
+     *
+     * In the multigrid preconditioner interfaces, the apply() method is used for
+     * the pre-smoothing operation because the previous content in the solution
+     * vector needs to be overwritten for a new incoming residual. On the other
+     * hand, all subsequent operations need to smooth the content already present
+     * in the vector @p u given the right hand side, which is done by smooth().
      */
     virtual void apply (const unsigned int level,
                         VectorType         &u,
@@ -352,7 +363,18 @@ public:
 
   /**
    * The apply variant of smoothing, setting the vector u to zero before
-   * calling the smooth function.
+   * calling the smooth function. This function is equivalent to the
+   * following code
+   * @code
+   * u = 0;
+   * smooth(level, u, rhs);
+   * @endcode
+   *
+   * In the multigrid preconditioner interfaces, the apply() method is used for
+   * the pre-smoothing operation because the previous content in the solution
+   * vector needs to be overwritten for a new incoming residual. On the other
+   * hand, all subsequent operations need to smooth the content already present
+   * in the vector @p u given the right hand side, which is done by smooth().
    */
   virtual void apply (const unsigned int level,
                       VectorType         &u,
@@ -493,7 +515,18 @@ public:
 
   /**
    * The apply variant of smoothing, setting the vector u to zero before
-   * calling the smooth function.
+   * calling the smooth function. This function is equivalent to the
+   * following code
+   * @code
+   * u = 0;
+   * smooth(level, u, rhs);
+   * @endcode
+   *
+   * In the multigrid preconditioner interfaces, the apply() method is used for
+   * the pre-smoothing operation because the previous content in the solution
+   * vector needs to be overwritten for a new incoming residual. On the other
+   * hand, all subsequent operations need to smooth the content already present
+   * in the vector @p u given the right hand side, which is done by smooth().
    */
   virtual void apply (const unsigned int level,
                       VectorType         &u,
index e79f7e0b20a0c41b84bac5d0a2d709e00204c25a..8f072ff8cd44ae1d79a65b874297aa9039b01052 100644 (file)
@@ -35,13 +35,9 @@ DEAL_II_NAMESPACE_OPEN
 /*@{*/
 
 /**
- * Implementation of the multigrid method.
- *
- * @warning multigrid on locally refined meshes only works with
- * <b>discontinuous finite elements</b> right now. It is not clear, whether
- * the paradigm of local smoothing we use is applicable to continuous elements
- * with hanging nodes; in fact, most people you meet on conferences seem to
- * deny this.
+ * Implementation of the multigrid method. The implementation supports both
+ * continuous and discontinuous elements and follows the procedure described in
+ * the @ref mg_paper "multigrid paper by Janssen and Kanschat".
  *
  * The function which starts a multigrid cycle on the finest level is cycle().
  * Depending on the cycle type chosen with the constructor (see enum Cycle),

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.