]> https://gitweb.dealii.org/ - dealii.git/commitdiff
new PreconditionRichardson
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 10 Mar 2005 17:30:34 +0000 (17:30 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 10 Mar 2005 17:30:34 +0000 (17:30 +0000)
git-svn-id: https://svn.dealii.org/trunk@10079 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.html
deal.II/lac/include/lac/precondition.h
tests/lac/solver.cc

index 70c631a8b5364f48326e6ae22846ae289ff9ebc1..a98486c5b1c519baf2b71a646c77999da246b645 100644 (file)
@@ -42,6 +42,13 @@ inconvenience this causes.
 
 
 <ol>
+  <li> <p>Changed: the file <tt>multigrid/mg_dof_tools.h</tt> was
+  renamed to <tt>multigrid/mg_tools.h</tt> in order to match the name
+  of the class contained.
+  <br>
+  (GK 2005/03/09)
+  </p>
+
   <li> <p>Changed: <code class="class">DoFTools</code>::<code
   class="member">make_flux_sparsity_pattern</code>, <code
   class="class">MGTools</code>::<code
@@ -157,39 +164,45 @@ inconvenience this causes.
 
 <ol>
   <li> <p>
-       Fixed: The <code>BlockSparseMatrix</code> class had no local
-       typedef <code>value_type</code> like all other classes, which
+       New: The <code class="class">PreconditionRichardson</code>
+       implements a Richardson preconditioner.
+
+       <br> 
+       (GK, 2005/03/10)
+       </p>
+
+  <li> <p> Fixed: The <code class="class">BlockSparseMatrix</code>
+       class had no local typedef <code
+       class="member">value_type</code> like all other classes, which
        made it a little awkward to use in some places. This has now
        been fixed.
+
        <br> 
        (WB, 2005/03/03)
        </p>
 
-  <li> <p>
-       Fixed: The <code>PETScWrappers::MatrixBase</code> class
-       documented that adding or setting a value that hasn't been in
-       the sparsity pattern before will lead to an exception being
-       thrown. This is of course wrong: PETSc allocates matrix entries
-       dynamically, as needed. The documentation is now fixed.
-       <br> 
-       (WB, 2005/03/03)
+  <li> <p> Fixed: The <code class="class">PETScWrappers</code>::<code
+       class="member">MatrixBase</code> class documented that adding
+       or setting a value that hasn't been in the sparsity pattern
+       before will lead to an exception being thrown. This is of
+       course wrong: PETSc allocates matrix entries dynamically, as
+       needed. The documentation is now fixed.
+       <br> (WB, 2005/03/03)
        </p>
 
-  <li> <p>
-       New: The <code>SparseMatrix</code> iterators had no <code>operator
-       &gt;</code>, only an <code>operator &lt;</code>. The missing operator
-       is now implemented. The same holds for the <code>FullMatrix</code>
-       class.
-       <br> 
-       (WB, 2005/03/01)
+  <li> <p> New: The <code class="class">SparseMatrix</code> iterators
+       had no <code class="member">operator &gt;</code>, only an <code
+       class="member">operator &lt;</code>. The missing operator is
+       now implemented. The same holds for the <code
+       class="class">FullMatrix</code> class.
+       <br> (WB, 2005/03/01)
        </p>
 
-  <li> <p>
-       Fixed: The <code>SparseMatrix</code> iterators could not be compared
-       using <code>operator &lt;</code>: the compiler complained that a
-       private member was accessed. This is now fixed.
-       <br> 
-       (WB, 2005/03/01)
+  <li> <p> Fixed: The <code class="class">SparseMatrix</code>
+       iterators could not be compared using <code
+       class="member">operator &lt;</code>: the compiler complained
+       that a private member was accessed. This is now fixed.
+       <br> (WB, 2005/03/01)
        </p>
 </ol>
 
index e71047fae8f7818b77c6b75f58c401e0bd2a920c..e4ee2156fa6f7f98a65b2027817bc4fd033485bc 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
+//    Copyright (C) 1998 - 2005 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -72,6 +72,97 @@ class PreconditionIdentity : public Subscriptor
 
 
 
+/**
+ * Preconditioning with Richardson's method. This preconditioner just
+ * scales the vector with a constant relaxation factor provided by the
+ * #AdditionalData object.
+ *
+ * In Krylov-space methods, this preconditioner should not have any
+ * effect. Using SolverRichardson, the two relaxation parameters will
+ * be just multiplied. Still, this class is usefull in multigrid
+ * smoother objects (MGSmootherRelaxation).
+ *
+ * @author Guido Kanschat, 2005
+ */
+class PreconditionRichardson : public Subscriptor
+{
+  public:
+                                    /**
+                                     * Parameters for Richardson
+                                     * preconditioner.
+                                     */
+    class AdditionalData
+    {
+      public:
+                                        /**
+                                         * Constructor. Block size
+                                         * must be given since there
+                                         * is no reasonable default
+                                         * parameter.
+                                         */
+       AdditionalData (const double relaxation = 1.);
+
+                                        /**
+                                         * Relaxation parameter.
+                                         */
+       double relaxation;
+    };
+
+                                    /**
+                                     * Constructor with optional
+                                     * setting of the relaxation
+                                     * parameter
+                                     */
+    PreconditionRichardson(const AdditionalData = AdditionalData());
+
+                                    /**
+                                     * Change the relaxaton parameter
+                                     * in a way consistent with other
+                                     * preconditioners.
+                                     */
+    void initialize (const AdditionalData parameters);
+
+                                    /**
+                                     * Apply preconditioner.
+                                     */
+    template<class VECTOR>
+    void vmult (VECTOR&, const VECTOR&) const;
+
+                                    /**
+                                     * Apply transpose
+                                     * preconditioner. Since this is
+                                     * the identity, this function is
+                                     * the same as
+                                     * vmult().
+                                     */
+    template<class VECTOR>
+    void Tvmult (VECTOR&, const VECTOR&) const;
+                                    /**
+                                     * Apply preconditioner, adding to the previous value.
+                                     */
+    template<class VECTOR>
+    void vmult_add (VECTOR&, const VECTOR&) const;
+
+                                    /**
+                                     * Apply transpose
+                                     * preconditioner, adding. Since this is
+                                     * the identity, this function is
+                                     * the same as
+                                     * vmult_add().
+                                     */
+    template<class VECTOR>
+    void Tvmult_add (VECTOR&, const VECTOR&) const;
+    
+  private:
+                                    /**
+                                     * The relaxation parameter
+                                     * multiplied with the vectors.
+                                     */
+    double relaxation;
+};
+
+
+
 /**
  * Preconditioner using a matrix-builtin function.
  * This class forms a preconditioner suitable for the LAC solver
@@ -641,6 +732,66 @@ PreconditionIdentity::Tvmult_add (VECTOR &dst, const VECTOR &src) const
 
 //----------------------------------------------------------------------//
 
+inline
+PreconditionRichardson::AdditionalData::AdditionalData (
+  const double relaxation)
+               :
+               relaxation(relaxation)
+{}
+
+
+inline
+PreconditionRichardson::PreconditionRichardson (
+  const PreconditionRichardson::AdditionalData parameters)
+               :
+               relaxation(parameters.relaxation)
+{}
+
+
+
+inline void
+PreconditionRichardson::initialize (
+  const PreconditionRichardson::AdditionalData parameters)
+{
+  relaxation = parameters.relaxation;
+}
+
+
+
+template<class VECTOR>
+inline void
+PreconditionRichardson::vmult (VECTOR &dst, const VECTOR &src) const
+{
+  dst.equ(relaxation,src);
+}
+
+
+
+template<class VECTOR>
+inline void
+PreconditionRichardson::Tvmult (VECTOR &dst, const VECTOR &src) const
+{
+  dst.equ(relaxation,src);
+}
+
+template<class VECTOR>
+inline void
+PreconditionRichardson::vmult_add (VECTOR &dst, const VECTOR &src) const
+{
+  dst.add(relaxation,src);
+}
+
+
+
+template<class VECTOR>
+inline void
+PreconditionRichardson::Tvmult_add (VECTOR &dst, const VECTOR &src) const
+{
+  dst.add(relaxation,src);
+}
+
+//----------------------------------------------------------------------//
+
 template <class MATRIX>
 inline void
 PreconditionRelaxation<MATRIX>::initialize (const MATRIX &rA,
index 3c905fc123272ce1bcb2275c3ac116f56a84b187..9c41ed07e81af22f89fabe87a8f2e1ac78311c47 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$ 
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
+//    Copyright (C) 1998 - 2005 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -99,6 +99,7 @@ int main()
       testproblem.five_point(A);
 
       PreconditionIdentity prec_no;
+      PreconditionRichardson prec_richardson(0.6);
       PreconditionSOR<> prec_sor;
       prec_sor.initialize(A, 1.2);
       PreconditionSSOR<> prec_ssor;
@@ -157,12 +158,28 @@ int main()
          deallog.pop();
          
          deallog.push("no");
-         
+
+         rich.set_omega(1./A.diag_element(0)); 
+         check_solve(rich,A,u,f,prec_no);
          check_solve(cg,A,u,f,prec_no);
          check_solve(bicgstab,A,u,f,prec_no);
          check_solve(gmres,A,u,f,prec_no);
          check_solve(gmresright,A,u,f,prec_no);
          check_solve(qmrs,A,u,f,prec_no);
+         rich.set_omega(1.);
+         
+         deallog.pop();
+         
+         deallog.push("rich");
+         
+         rich.set_omega(1./A.diag_element(0)); 
+         check_solve(rich,A,u,f,prec_richardson);
+         check_solve(cg,A,u,f,prec_richardson);
+         check_solve(bicgstab,A,u,f,prec_richardson);
+         check_solve(gmres,A,u,f,prec_richardson);
+         check_solve(gmresright,A,u,f,prec_richardson);
+         check_solve(qmrs,A,u,f,prec_richardson);
+         rich.set_omega(1.);
          
          deallog.pop();
          

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.