]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Vanka works
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 5 Jul 1999 13:48:18 +0000 (13:48 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 5 Jul 1999 13:48:18 +0000 (13:48 +0000)
git-svn-id: https://svn.dealii.org/trunk@1536 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/sparse_matrix.2.templates
deal.II/lac/include/lac/sparse_vanka.h [moved from deal.II/lac/include/lac/sparsevanka.h with 66% similarity]
deal.II/lac/include/lac/sparse_vanka.templates.h [moved from deal.II/lac/include/lac/sparsevanka.templates.h with 80% similarity]
deal.II/lac/source/sparse_matrix.double.cc
deal.II/lac/source/sparse_matrix.float.cc
deal.II/lac/source/sparse_matrix.long_double.cc

index 177b86395dce1f44eda81d88fa1b0b2530b6de0d..1be0e1294e84cbe028adf8013a213c300a07db9d 100644 (file)
@@ -43,5 +43,6 @@ template void SparseMatrix<TYPEMAT>::SOR_step (Vector<TYPE2> &, const Vector<TYP
 template void SparseMatrix<TYPEMAT>::TSOR_step (Vector<TYPE2> &, const Vector<TYPE2> &, TYPEMAT) const;
 template void SparseMatrix<TYPEMAT>::SSOR_step (Vector<TYPE2> &, const Vector<TYPE2> &, TYPEMAT) const;
 
-template void SparseVanka<TYPEMAT>::apply(Vector<TYPE2>& dst) const;
+template void SparseVanka<TYPEMAT>::forward(Vector<TYPE2>& dst, const Vector<TYPE2>& dst) const;
+template void SparseVanka<TYPEMAT>::backward(Vector<TYPE2>& dst, const Vector<TYPE2>& dst) const;
 
similarity index 66%
rename from deal.II/lac/include/lac/sparsevanka.h
rename to deal.II/lac/include/lac/sparse_vanka.h
index 0a0664287a42ce98cbc96007bd64200c1b1a29af..9b3537280eeb66e8de989c73b3ddd1ee1cf30042 100644 (file)
@@ -1,8 +1,8 @@
 // $Id$
 // Copyright Guido Kanschat, 1999
 
-#ifndef __lac_sparsematrix_H
-#define __lac_sparsematrix_H
+#ifndef __lac_sparse_vanka_H
+#define __lac_sparse_vanka_H
 
 #include <base/smartpointer.h>
 #include <lac/forward-declarations.h>
@@ -29,8 +29,9 @@
  *
  * This local system is solved and the values are updated into the
  * destination vector.
- * @author Guido Kanschat
- */
+ *
+ * Remark: the Vanka method is a non-symmetric preconditioning method.
+ * @author Guido Kanschat */
 template<typename number>
 class SparseVanka
 {
@@ -48,19 +49,36 @@ class SparseVanka
     SparseVanka(const SparseMatrix<number>& M,
                const vector<int>& indices);
                                     /**
-                                     * Do the preconditioning.
+                                     * Do the preconditioning. This
+                                     * function contains a dispatch
+                                     * mechanism to use the
+                                     * multiplicative version by
+                                     * default and the additive version
+                                     * if requested by #set_additive#.
                                      */
     template<typename number2>
     void operator() (Vector<number2>& dst,
                     const Vector<number2>& src) const;
 
                                     /**
-                                     * In-place application of the
-                                     * method.
+                                     * Application of the Vanka operator.
+                                     * This function takes two vector
+                                     * arguments, the residual in #src#
+                                     * and the resulting update vector
+                                     * in #dst#.
                                      */
     template<typename number2>
-    void apply(Vector<number2>& dst) const;
-    
+    void forward(Vector<number2>& dst, const Vector<number2>& src) const;
+                                    /**
+                                     * Application of the transpose
+                                     * Vanka operator.
+                                     * This function takes two vector
+                                     * arguments, the residual in #src#
+                                     * and the resulting update vector
+                                     * in #dst#.
+                                     */
+    template<typename number2>
+    void backward(Vector<number2>& dst, const Vector<number2>& src) const;
   private:
                                     /**
                                      * Pointer to the matrix.
@@ -81,8 +99,8 @@ void
 SparseVanka<number>::operator() (Vector<number2>& dst,
                                 const Vector<number2>& src) const
 {
-  dst = src;
-  apply(dst);
+  dst = 0.;
+  forward(dst, src);
 }
 
 
similarity index 80%
rename from deal.II/lac/include/lac/sparsevanka.templates.h
rename to deal.II/lac/include/lac/sparse_vanka.templates.h
index 19e87b7d60cec3efd77ea2e35af025ce53492513..be64dbd95560615144caaad1cf477271dccb923f 100644 (file)
@@ -1,7 +1,7 @@
 // $Id$
 // Copyright Guido Kanschat, 1999
 
-#include <lac/sparsevanka.h>
+#include <lac/sparse_vanka.h>
 #include <lac/fullmatrix.h>
 
 #include <map>
@@ -16,7 +16,8 @@ SparseVanka<number>::SparseVanka(const SparseMatrix<number>& M,
 template<typename number>
 template<typename number2>
 void
-SparseVanka<number>::apply(Vector<number2>& dst) const
+SparseVanka<number>::forward(Vector<number2>& dst,
+                          const Vector<number2>& src) const
 {
   for (unsigned int global_i=0; global_i<indices.size() ; ++global_i)
     {
@@ -36,6 +37,10 @@ SparseVanka<number>::apply(Vector<number2>& dst) const
        local_index.insert(pair<unsigned int, unsigned int>
                           (structure.column_number(row, i), i));
 
+//       for (map<unsigned int, unsigned int>::iterator is=local_index.begin();
+//        is!=local_index.end();++is)
+//     cerr << "map " << is->first << '\t' << is->second << endl;
+      
                                       // Build local matrix
 
       for (map<unsigned int, unsigned int>::iterator is=local_index.begin();
@@ -45,7 +50,7 @@ SparseVanka<number>::apply(Vector<number2>& dst) const
          unsigned int i = is->second;
          unsigned int n = structure.row_length(irow);
          
-         b(i) = dst(irow);
+         b(i) = src(irow);
          
          for (unsigned int j=0;j<n;++j)
            {
@@ -54,9 +59,9 @@ SparseVanka<number>::apply(Vector<number2>& dst) const
                = local_index.find(col);
              if (js == local_index.end())
                {
-                 b(i) -= matrix->raw_entry(irow,col) * dst(col);
+                 b(i) -= matrix->raw_entry(irow,j) * dst(col);
                } else {
-                 A(i,j) = matrix->raw_entry(irow,col);
+                 A(i,js->second) = matrix->raw_entry(irow,j);
                }
            }
        }
@@ -71,7 +76,7 @@ SparseVanka<number>::apply(Vector<number2>& dst) const
          unsigned int irow = is->first;
          unsigned int i = is->second;
          dst(irow) = x(i);
-       }     
+       }
     }
 }
 
index d677b3275fe48bc36442e086abfe4f23b0fc293b..552653314c4abef52bf91aa055e4269ac75722ca 100644 (file)
@@ -11,7 +11,7 @@
 
 #include <cmath>
 #include <lac/sparsematrix.templates.h>
-#include <lac/sparsevanka.templates.h>
+#include <lac/sparse_vanka.templates.h>
 
 #define TYPEMAT double
 
index 32b2e0324e50d96bd413bc5618b035891911ef90..996aa409d6e4b8f31b3b35b97c55b546ac63ab83 100644 (file)
@@ -11,7 +11,7 @@
 
 #include <cmath>
 #include <lac/sparsematrix.templates.h>
-#include <lac/sparsevanka.templates.h>
+#include <lac/sparse_vanka.templates.h>
 
 #define TYPEMAT float
 
index 9afdba45309ae5a6d9334d1146488b7268734938..4a1b619177889d483a0f76afc108e23cb62a7655 100644 (file)
@@ -11,7 +11,7 @@
 
 #include <cmath>
 #include <lac/sparsematrix.templates.h>
-#include <lac/sparsevanka.templates.h>
+#include <lac/sparse_vanka.templates.h>
 
 #define TYPEMAT long double
 

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.