]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
unifying preconditioners
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 6 Jan 2003 14:55:47 +0000 (14:55 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 6 Jan 2003 14:55:47 +0000 (14:55 +0000)
git-svn-id: https://svn.dealii.org/trunk@6868 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/precondition.h
deal.II/lac/include/lac/precondition_block.h
deal.II/lac/include/lac/precondition_block.templates.h

index a8cc49fe601c397aa3b0ecd04886b89063188a10..62e6652285944f1de7fa9b9262a0b972f2f7f35e 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002 by the deal.II authors
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -167,6 +167,23 @@ template<class MATRIX = SparseMatrix<double> >
 class PreconditionRelaxation : public Subscriptor
 {
   public:
+                                    /**
+                                     * Class for parameters.
+                                     */
+    class AdditionalData
+    {
+      public:
+                                        /**
+                                         * Constructor.
+                                         */
+       AdditionalData (const double relaxation = 1.);
+
+                                        /**
+                                         * Relaxation parameter.
+                                         */
+       double relaxation;      
+    };
+    
                                     /**
                                      * Initialize matrix and
                                      * relaxation parameter. The
@@ -177,7 +194,8 @@ class PreconditionRelaxation : public Subscriptor
                                      * than 2 for numerical
                                      * reasons. It defaults to 1.
                                      */
-    void initialize (const MATRIX& A, const double omega = 1.);
+    void initialize (const MATRIX& A,
+                    AdditionalData parameters = AdditionalData());
     
   protected:
                                     /**
@@ -188,7 +206,7 @@ class PreconditionRelaxation : public Subscriptor
                                     /**
                                      * Relaxation parameter.
                                      */
-    double omega;
+    double relaxation;
 };
 
 
@@ -531,10 +549,10 @@ PreconditionIdentity::Tvmult_add (VECTOR& dst, const VECTOR& src) const
 template <class MATRIX>
 inline void
 PreconditionRelaxation<MATRIX>::initialize (const MATRIX &rA,
-                                           const double  o)
+                                           AdditionalData parameters)
 {
   A = &rA;
-  omega = o;
+  relaxation = parameters.relaxation;
 }
 
 //----------------------------------------------------------------------//
@@ -545,7 +563,7 @@ inline void
 PreconditionJacobi<MATRIX>::vmult (VECTOR& dst, const VECTOR& src) const
 {
   Assert (this->A!=0, ExcNotInitialized());
-  this->A->precondition_Jacobi (dst, src, this->omega);
+  this->A->precondition_Jacobi (dst, src, this->relaxation);
 }
 
 
@@ -556,7 +574,7 @@ inline void
 PreconditionJacobi<MATRIX>::Tvmult (VECTOR& dst, const VECTOR& src) const
 {
   Assert (this->A!=0, ExcNotInitialized());
-  this->A->precondition_Jacobi (dst, src, this->omega);
+  this->A->precondition_Jacobi (dst, src, this->relaxation);
 }
 
 
@@ -568,7 +586,7 @@ inline void
 PreconditionSOR<MATRIX>::vmult (VECTOR& dst, const VECTOR& src) const
 {
   Assert (this->A!=0, ExcNotInitialized());
-  this->A->precondition_SOR (dst, src, this->omega);
+  this->A->precondition_SOR (dst, src, this->relaxation);
 }
 
 
@@ -579,7 +597,7 @@ inline void
 PreconditionSOR<MATRIX>::Tvmult (VECTOR& dst, const VECTOR& src) const
 {
   Assert (this->A!=0, ExcNotInitialized());
-  this->A->precondition_TSOR (dst, src, this->omega);
+  this->A->precondition_TSOR (dst, src, this->relaxation);
 }
 
 
@@ -591,7 +609,7 @@ inline void
 PreconditionSSOR<MATRIX>::vmult (VECTOR& dst, const VECTOR& src) const
 {
   Assert (this->A!=0, ExcNotInitialized());
-  this->A->precondition_SSOR (dst, src, this->omega);
+  this->A->precondition_SSOR (dst, src, this->relaxation);
 }
 
 
@@ -602,7 +620,7 @@ inline void
 PreconditionSSOR<MATRIX>::Tvmult (VECTOR& dst, const VECTOR& src) const
 {
   Assert (this->A!=0, ExcNotInitialized());
-  this->A->precondition_SSOR (dst, src, this->omega);
+  this->A->precondition_SSOR (dst, src, this->relaxation);
 }
 
 
@@ -626,6 +644,18 @@ PreconditionUseMatrix<MATRIX,VECTOR>::vmult (VECTOR& dst,
   (matrix.*precondition)(dst, src);
 }
 
+//----------------------------------------------------------------------//
+
+template<class MATRIX>
+inline
+PreconditionRelaxation<MATRIX>::AdditionalData::
+AdditionalData (const double relaxation)
+               :
+               relaxation (relaxation)
+{}
+
+
+
 //////////////////////////////////////////////////////////////////////
 
 template<class SOLVER, class MATRIX, class PRECONDITION>
index 37dab620a85fb155a1325f6b81656b8144847178..4d2f8d2ffca41fbe3e207f895fa6912283a4813a 100644 (file)
@@ -85,6 +85,33 @@ class PreconditionBlock : public virtual Subscriptor
     typedef typename MATRIX::value_type number;
     
   public:
+                                    /**
+                                     * Parameters for block preconditioners.
+                                     */
+    class AdditionalData
+    {
+      public:
+                                        /**
+                                         * Constructor. Block size
+                                         * must be given since there
+                                         * is no reasonable default
+                                         * parameter.
+                                         */
+       AdditionalData (const unsigned int block_size,
+                       const double relaxation = 1.);
+
+                                        /**
+                                         * Relaxation parameter.
+                                         */
+       double relaxation;
+       
+                                        /**
+                                         * Block size.
+                                         */
+       unsigned int block_size;
+    };
+    
+    
                                     /**
                                      * Constructor.
                                      */
@@ -108,9 +135,7 @@ class PreconditionBlock : public virtual Subscriptor
                                      * parameter for derived classes
                                      * may be provided.
                                      */
-    void initialize (const MATRIX& A,
-                    const unsigned int block_size,
-                    const double relaxation = 1.);
+    void initialize (const MATRIX& A, AdditionalData parameters);
 
                                     /**
                                      * Deletes the inverse diagonal
@@ -556,6 +581,17 @@ class PreconditionBlockSSOR : public virtual Subscriptor,
 
 //----------------------------------------------------------------------//
 
+template<class MATRIX, typename inverse_type>
+inline
+PreconditionBlock<MATRIX, inverse_type>::AdditionalData::
+AdditionalData (const unsigned int block_size,
+               const double relaxation)
+               :
+               relaxation (relaxation),
+               block_size(block_size)
+{}
+
+
 template<class MATRIX, typename inverse_type>
 inline bool
 PreconditionBlock<MATRIX, inverse_type>::empty () const
index bac456a0511e32a59c0eba72ef06a46f4ad2f6b7..790cc4f2ff19b4aaede73caeab9abcfd845ae821 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002 by the deal.II authors
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -57,17 +57,19 @@ void PreconditionBlock<MATRIX,inverse_type>::clear ()
 
 
 template <class MATRIX, typename inverse_type>
-void PreconditionBlock<MATRIX,inverse_type>::initialize (const MATRIX &M,
-                                                        const unsigned int bsize,
-                                                        const double relax)
+void PreconditionBlock<MATRIX,inverse_type>::initialize (
+  const MATRIX &M,
+  AdditionalData parameters)
 {
+  const unsigned int bsize = parameters.block_size;
+  
   clear();
   Assert (M.m() == M.n(), ExcMatrixNotSquare());
   A = &M;
   Assert (bsize>0, ExcIndexRange(bsize, 1, M.m()));
   Assert (A->m()%bsize==0, ExcWrongBlockSize(bsize, A->m()));
   blocksize=bsize;
-  relaxation = relax;
+  relaxation = parameters.relaxation;
 }
 
 

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.