]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Make the print_vectors function usable when MG_DEBUG is defined. Implementation of...
authorhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 3 Feb 2003 10:34:30 +0000 (10:34 +0000)
committerhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 3 Feb 2003 10:34:30 +0000 (10:34 +0000)
git-svn-id: https://svn.dealii.org/trunk@7003 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/multigrid/multigrid.h
deal.II/deal.II/include/multigrid/multigrid.templates.h
deal.II/deal.II/source/multigrid/multigrid.all_dimensions.cc

index 06368d3a47c6b776eeb7dc8ba36426833710e382..fc924440cab35a162ba5755e44bd084cc3ecbe4d 100644 (file)
@@ -59,7 +59,7 @@
  * returning the solution to the system @p{Ax=b} on the coarsest level
  * in @p{x}.
  * 
- * @author Guido Kanschat, 1999, 2001, 2002, Ralf Hartmann 2002
+ * @author Guido Kanschat, 1999, 2001, 2002, Ralf Hartmann 2002, 2003.
  */
 template <class VECTOR>
 class Multigrid : public Subscriptor
@@ -137,6 +137,13 @@ class Multigrid : public Subscriptor
                                      */
     void set_edge_matrices (const MGMatrixBase<VECTOR>& edge_down,
                            const MGMatrixBase<VECTOR>& edge_up);
+
+#ifdef MG_DEBUG
+    void print_vector (const unsigned int level,
+                      const VECTOR &v,
+                      const char *name) const;
+#endif
+    
   private:
     
                                     /**
@@ -154,6 +161,7 @@ class Multigrid : public Subscriptor
                                      * this level.
                                      */
     void level_mgstep (const unsigned int level);
+
                                     /**
                                      * Level for coarse grid solution.
                                      */
@@ -220,6 +228,15 @@ class Multigrid : public Subscriptor
                                      */
     SmartPointer<const MGMatrixBase<VECTOR> > edge_up;
 
+                                    /**
+                                     * Pointer to the MGDoFHandler
+                                     * given to the constructor. Only
+                                     * needed for @p{MG_DEBUG}
+                                     * defined.
+                                     */
+#ifdef MG_DEBUG
+    SmartPointer<const MGDoFHandler<deal_II_dimension> > mg_dof_handler;
+#endif
                                     /**
                                      * Exception.
                                      */
@@ -345,6 +362,9 @@ Multigrid<VECTOR>::Multigrid (const MGDoFHandler<dim>& mg_dof_handler,
                post_smooth(&post_smooth),
                edge_down(0),
                edge_up(0)
+#ifdef MG_DEBUG  
+                , mg_dof_handler(&mg_dof_handler)
+#endif
 {}
 
 
index 60e00fddcf2a21d15db30f379039730f60d31ecf..e235a3a135f6a8e26ae8066a84c3f976418ab644 100644 (file)
 //----------------------------  multigrid.templates.h  ---------------------------
 #ifndef __deal2__multigrid_templates_h
 #define __deal2__multigrid_templates_h
+#include <multigrid/multigrid.h>
+
+#ifdef MG_DEBUG
+#include <fe/fe.h>
+#include <dofs/dof_accessor.h>
+#include <numerics/data_out.h>
+#include <multigrid/mg_dof_accessor.h>
+#include <multigrid/mg_dof_handler.h>
+#include <fstream>
+#endif
+
 
 
-#include <multigrid/multigrid.h>
 
 template <class VECTOR>
 void
@@ -50,7 +60,8 @@ Multigrid<VECTOR>::level_mgstep(const unsigned int level)
       return;
     }
 
-                          // smoothing of the residual by modifying s
+                                  // smoothing of the residual by
+                                  // modifying s
   pre_smooth->smooth(level, solution[level], defect[level]);
 
 #ifdef MG_DEBUG
@@ -61,11 +72,11 @@ Multigrid<VECTOR>::level_mgstep(const unsigned int level)
                                   // t = A*solution[level]
   matrix->vmult(level, t[level], solution[level]);
   
-                                  // make t rhs of lower level
-                                  // The non-refined parts of the
-                                  // coarse-level defect already contain
-                                  // the global defect, the refined parts
-                                  // its restriction.
+                                  // make t rhs of lower level The
+                                  // non-refined parts of the
+                                  // coarse-level defect already
+                                  // contain the global defect, the
+                                  // refined parts its restriction.
   for (unsigned int l = level;l>minlevel;--l)
     {
       t[l-1] = 0.;
@@ -92,7 +103,7 @@ Multigrid<VECTOR>::level_mgstep(const unsigned int level)
 
 #ifdef MG_DEBUG
   sprintf(name, "MG%d-cgc",level);
-  print_vector(level, t, name);
+  print_vector(level, t[level], name);
 #endif
 
   solution[level] += t[level];
@@ -137,54 +148,5 @@ Multigrid<VECTOR>::vcycle()
 
 
 
-// template <int dim>
-// void
-// Multigrid<dim>::print_vector (const unsigned int /*level*/,
-//                           const Vector<double>& /*v*/,
-//                           const char* /*name*/) const
-// {
-//    if (level!=maxlevel)
-//      return;
-  
-//    const DoFHandler<dim>* dof = mg_dof_handler;
-  
-//    Vector<double> out_vector(dof->n_dofs());
-
-//    out_vector = -10000;
-  
-//    const unsigned int dofs_per_cell = mg_dof_handler->get_fe().dofs_per_cell;
-
-//    std::vector<unsigned int> global_dof_indices (dofs_per_cell);
-//    std::vector<unsigned int> level_dof_indices (dofs_per_cell);
-
-//    typename DoFHandler<dim>::active_cell_iterator
-//      global_cell = dof->begin_active(level);
-//    typename MGDoFHandler<dim>::active_cell_iterator
-//      level_cell = mg_dof_handler->begin_active(level);
-//    const typename MGDoFHandler<dim>::cell_iterator
-//      endc = mg_dof_handler->end(level);
-
-//                                // traverse all cells and copy the
-//                                // data appropriately to the output
-//                                // vector
-//    for (; level_cell != endc; ++level_cell, ++global_cell)
-//      {
-//        global_cell->get_dof_indices (global_dof_indices);
-//        level_cell->get_mg_dof_indices(level_dof_indices);
-
-//                                    // copy level-wise data to
-//                                    // global vector
-//        for (unsigned int i=0; i<dofs_per_cell; ++i)
-//     out_vector(global_dof_indices[i])
-//       = v(level_dof_indices[i]);
-//      }
-
-//    std::ofstream out_file(name);
-//    DataOut<dim> out;
-//    out.attach_dof_handler(*dof);
-//    out.add_data_vector(out_vector, "v");
-//    out.build_patches(5);
-//    out.write_gnuplot(out_file);
-//}
 
 #endif
index f734dfc6d56b900a8e15cc8ba5bb1394d03df08f..b634c62f5eed85b467771b5c9828abafb48abb19 100644 (file)
 #include <multigrid/multigrid.templates.h>
 
 
+
+
+
+#ifdef MG_DEBUG
+template <>
+void
+Multigrid<Vector<double> >::print_vector (const unsigned int level,
+                                         const Vector<double> &v,
+                                         const char *name) const
+{
+   if (level!=maxlevel)
+     return;
+   const unsigned int dim=deal_II_dimension;
+   
+   const DoFHandler<dim> *dof = mg_dof_handler;
+  
+   Vector<double> out_vector(dof->n_dofs());
+
+   out_vector = -10000;
+  
+   const unsigned int dofs_per_cell = mg_dof_handler->get_fe().dofs_per_cell;
+
+   std::vector<unsigned int> global_dof_indices (dofs_per_cell);
+   std::vector<unsigned int> level_dof_indices (dofs_per_cell);
+
+   DoFHandler<dim>::active_cell_iterator
+     global_cell = dof->begin_active(level);
+   MGDoFHandler<dim>::active_cell_iterator
+     level_cell = mg_dof_handler->begin_active(level);
+   const MGDoFHandler<dim>::cell_iterator
+     endc = mg_dof_handler->end(level);
+
+                                  // traverse all cells and copy the
+                                  // data appropriately to the output
+                                  // vector
+   for (; level_cell != endc; ++level_cell, ++global_cell)
+     {
+       global_cell->get_dof_indices (global_dof_indices);
+       level_cell->get_mg_dof_indices(level_dof_indices);
+
+                                      // copy level-wise data to
+                                      // global vector
+       for (unsigned int i=0; i<dofs_per_cell; ++i)
+       out_vector(global_dof_indices[i])
+         = v(level_dof_indices[i]);
+     }
+
+   std::ofstream out_file(name);
+   DataOut<dim> out;
+   out.attach_dof_handler(*dof);
+   out.add_data_vector(out_vector, "v");
+   out.build_patches(5);
+   out.write_gnuplot(out_file);
+}
+
+
+
+template <class VECTOR>
+void
+Multigrid<VECTOR>::print_vector (const unsigned int level,
+                                const VECTOR &v,
+                                const char *name) const
+{
+  Assert(false, ExcNotImplemented());
+}
+#endif
+
+
+
+
 template <typename number>
 MGTransferPrebuilt<number>::~MGTransferPrebuilt () 
 {}

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.