]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
add coarse grid solver based on SVD
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 3 Apr 2012 18:45:11 +0000 (18:45 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 3 Apr 2012 18:45:11 +0000 (18:45 +0000)
git-svn-id: https://svn.dealii.org/trunk@25377 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/multigrid/mg_coarse.h

index 7a2dba40359f8014570603df4a44c105162b5942..afb28b8dfcbafb4c70cc48afe8e751a1968b7cf5 100644 (file)
@@ -134,25 +134,58 @@ class MGCoarseGridHouseholder : public MGCoarseGridBase<VECTOR>
                                       * Initialize for a new matrix.
                                       */
     void initialize (const FullMatrix<number>& A);
+    
+    void operator() (const unsigned int   level,
+                     VECTOR       &dst,
+                     const VECTOR &src) const;
 
+  private:
                                      /**
-                                      * Release matrix pointer.
+                                      * Matrix for QR-factorization.
                                       */
-    void clear ();
+    Householder<number> householder;
+};
 
+/**
+ * Coarse grid solver using singular value decomposition of LAPACK matrices.
+ *
+ * Upon initialization, the singular value decomposition of the matrix is
+ * computed. then, the operator() uses 
+ *
+ * @author Guido Kanschat, 2003, 2012
+ */
+template<typename number = double, class VECTOR = Vector<number> >
+class MGCoarseGridSVD : public MGCoarseGridBase<VECTOR>
+{
+  public:
+                                     /**
+                                      * Constructor leaving an
+                                      * uninitialized object.
+                                      */
+    MGCoarseGridSVD ();
+    
                                      /**
-                                      * Solution operator, defined in
-                                      * the base class.
+                                      * Initialize for a new
+                                      * matrix. This resets the
+                                      * dimensions to the 
                                       */
+    void initialize (const FullMatrix<number>& A, const double threshold = 1.e-12);
+    
     void operator() (const unsigned int   level,
                      VECTOR       &dst,
                      const VECTOR &src) const;
-
+    
+                                    /**
+                                     * Write the singular values to #deallog.
+                                     */
+    void log () const;
+    
   private:
+    
                                      /**
-                                      * Matrix for QR-factorization.
+                                      * Matrix for singular value decomposition.
                                       */
-    Householder<number> householder;
+    LAPACKFullMatrix<number> matrix;
 };
 
 /*@}*/
@@ -274,21 +307,58 @@ MGCoarseGridHouseholder<number, VECTOR>::initialize(
 
 template<typename number, class VECTOR>
 void
-MGCoarseGridHouseholder<number, VECTOR>::clear()
+MGCoarseGridHouseholder<number, VECTOR>::operator() (
+  const unsigned int,
+  VECTOR       &dst,
+  const VECTOR &src) const
+{
+  householder.least_squares(dst, src);
+}
+
+//---------------------------------------------------------------------------
+
+template<typename number, class VECTOR>
+inline
+MGCoarseGridSVD<number, VECTOR>::MGCoarseGridSVD()
 {}
 
 
 
 template<typename number, class VECTOR>
 void
-MGCoarseGridHouseholder<number, VECTOR>::operator() (
+MGCoarseGridSVD<number, VECTOR>::initialize(
+  const FullMatrix<number>& A,
+  double threshold)
+{
+  matrix.reinit(A.n_rows(), A.n_cols());
+  matrix = A;
+  matrix.compute_inverse_svd(threshold);
+}
+
+
+template<typename number, class VECTOR>
+void
+MGCoarseGridSVD<number, VECTOR>::operator() (
   const unsigned int,
   VECTOR       &dst,
   const VECTOR &src) const
 {
-  householder.least_squares(dst, src);
+  matrix.vmult(dst, src);
 }
 
+
+template<typename number, class VECTOR>
+void
+MGCoarseGridSVD<number, VECTOR>::log() const
+{
+  const unsigned int n = std::min(matrix.n_rows(), matrix.n_cols());
+  
+  for (unsigned int i=0;i<n;++i)
+    deallog << ' ' << matrix.singular_value(i);
+  deallog << std::endl;
+}
+
+
 #endif // DOXYGEN
 
 DEAL_II_NAMESPACE_CLOSE

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.