]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove Multigrid::set_debug
authorDaniel Arndt <arndtd@ornl.gov>
Sun, 29 Mar 2020 13:37:39 +0000 (09:37 -0400)
committerDaniel Arndt <arndtd@ornl.gov>
Sun, 29 Mar 2020 13:41:36 +0000 (09:41 -0400)
include/deal.II/multigrid/multigrid.h
include/deal.II/multigrid/multigrid.templates.h

index 32f54b6a7f7b2d2a2803bdae173ed930cdda6a29..dec545aa68f044b32f307993a0e8605426141fad 100644 (file)
@@ -272,16 +272,6 @@ public:
    */
   void set_cycle(Cycle);
 
-  /**
-   * @deprecated Debug output will go away. Use signals instead.
-   *
-   * Set the debug level. Higher values will create more debugging output
-   * during the multigrid cycles.
-   */
-  DEAL_II_DEPRECATED
-  void
-  set_debug(const unsigned int);
-
   /**
    * Connect a function to mg::Signals::coarse_solve.
    */
@@ -439,11 +429,6 @@ private:
    */
   SmartPointer<const MGMatrixBase<VectorType>, Multigrid<VectorType>> edge_up;
 
-  /**
-   * Level for debug output. Defaults to zero and can be set by set_debug().
-   */
-  unsigned int debug;
-
   template <int dim, class OtherVectorType, class TRANSFER>
   friend class PreconditionMG;
 };
@@ -626,7 +611,6 @@ Multigrid<VectorType>::Multigrid(const MGMatrixBase<VectorType> &    matrix,
   , edge_in(nullptr, typeid(*this).name())
   , edge_down(nullptr, typeid(*this).name())
   , edge_up(nullptr, typeid(*this).name())
-  , debug(0)
 {
   if (max_level == numbers::invalid_unsigned_int)
     maxlevel = matrix.get_maxlevel();
index 108f3ab6547a475bc6e117cb4f38b76dc52a5fed..c21e1f5ff5f13a1807872e2839a1a77036927f46 100644 (file)
@@ -74,15 +74,6 @@ Multigrid<VectorType>::set_cycle(typename Multigrid<VectorType>::Cycle c)
 
 
 
-template <typename VectorType>
-void
-Multigrid<VectorType>::set_debug(const unsigned int d)
-{
-  debug = d;
-}
-
-
-
 template <typename VectorType>
 void
 Multigrid<VectorType>::set_edge_matrices(const MGMatrixBase<VectorType> &down,
@@ -110,50 +101,27 @@ template <typename VectorType>
 void
 Multigrid<VectorType>::level_v_step(const unsigned int level)
 {
-  if (debug > 0)
-    deallog << "V-cycle entering level " << level << std::endl;
-  if (debug > 2)
-    deallog << "V-cycle  Defect norm   " << defect[level].l2_norm()
-            << std::endl;
-
   if (level == minlevel)
     {
       this->signals.coarse_solve(true, level);
-      if (debug > 0)
-        deallog << "Coarse level           " << level << std::endl;
       (*coarse)(level, solution[level], defect[level]);
       this->signals.coarse_solve(false, level);
       return;
     }
 
   // smoothing of the residual
-  if (debug > 1)
-    deallog << "Smoothing on     level " << level << std::endl;
-
   this->signals.pre_smoother_step(true, level);
   pre_smooth->apply(level, solution[level], defect[level]);
   this->signals.pre_smoother_step(false, level);
 
-  if (debug > 2)
-    deallog << "Solution norm          " << solution[level].l2_norm()
-            << std::endl;
-
   // compute residual on level, which includes the (CG) edge matrix
-  if (debug > 1)
-    deallog << "Residual on      level " << level << std::endl;
   matrix->vmult(level, t[level], solution[level]);
   if (edge_out != nullptr)
     {
       edge_out->vmult_add(level, t[level], solution[level]);
-      if (debug > 2)
-        deallog << "Norm     t[" << level << "] " << t[level].l2_norm()
-                << std::endl;
     }
   t[level].sadd(-1.0, 1.0, defect[level]);
 
-  if (debug > 2)
-    deallog << "Residual norm          " << t[level].l2_norm() << std::endl;
-
   // Get the defect on the next coarser level as part of the (DG) edge matrix
   // and then the main part by the restriction of the transfer
   if (edge_down != nullptr)
@@ -174,8 +142,6 @@ Multigrid<VectorType>::level_v_step(const unsigned int level)
   transfer->prolongate(level, t[level], solution[level - 1]);
   this->signals.prolongation(false, level);
 
-  if (debug > 2)
-    deallog << "Prolongate norm        " << t[level].l2_norm() << std::endl;
   solution[level] += t[level];
 
   // get in contribution from edge matrices to the defect
@@ -190,23 +156,10 @@ Multigrid<VectorType>::level_v_step(const unsigned int level)
       defect[level] -= t[level];
     }
 
-  if (debug > 2)
-    deallog << "V-cycle  Defect norm   " << defect[level].l2_norm()
-            << std::endl;
-
   // post-smoothing
-  if (debug > 1)
-    deallog << "Smoothing on     level " << level << std::endl;
   this->signals.post_smoother_step(true, level);
   post_smooth->smooth(level, solution[level], defect[level]);
   this->signals.post_smoother_step(false, level);
-
-  if (debug > 2)
-    deallog << "Solution norm          " << solution[level].l2_norm()
-            << std::endl;
-
-  if (debug > 1)
-    deallog << "V-cycle leaving  level " << level << std::endl;
 }
 
 
@@ -231,48 +184,26 @@ Multigrid<VectorType>::level_step(const unsigned int level, Cycle cycle)
         Assert(false, ExcNotImplemented());
     }
 
-  if (debug > 0)
-    deallog << cychar << "-cycle entering level  " << level << std::endl;
-
   // Combine the defect from the initial copy_to_mg with the one that has come
   // from the finer level by the transfer
   defect2[level] += defect[level];
   defect[level] = typename VectorType::value_type(0.);
 
-  if (debug > 2)
-    deallog << cychar << "-cycle defect norm     " << defect2[level].l2_norm()
-            << std::endl;
-
   if (level == minlevel)
     {
-      if (debug > 0)
-        deallog << cychar << "-cycle coarse level    " << level << std::endl;
-
       (*coarse)(level, solution[level], defect2[level]);
       return;
     }
 
   // smoothing of the residual
-  if (debug > 1)
-    deallog << cychar << "-cycle smoothing level " << level << std::endl;
   pre_smooth->apply(level, solution[level], defect2[level]);
 
-  if (debug > 2)
-    deallog << cychar << "-cycle solution norm   " << solution[level].l2_norm()
-            << std::endl;
-
   // compute residual on level, which includes the (CG) edge matrix
-  if (debug > 1)
-    deallog << cychar << "-cycle residual level  " << level << std::endl;
   matrix->vmult(level, t[level], solution[level]);
   if (edge_out != nullptr)
     edge_out->vmult_add(level, t[level], solution[level]);
   t[level].sadd(-1.0, 1.0, defect2[level]);
 
-  if (debug > 2)
-    deallog << cychar << "-cycle residual norm   " << t[level].l2_norm()
-            << std::endl;
-
   // Get the defect on the next coarser level as part of the (DG) edge matrix
   // and then the main part by the restriction of the transfer
   if (edge_down != nullptr)
@@ -315,21 +246,8 @@ Multigrid<VectorType>::level_step(const unsigned int level, Cycle cycle)
       defect2[level] -= t[level];
     }
 
-  if (debug > 2)
-    deallog << cychar << "-cycle  Defect norm    " << defect2[level].l2_norm()
-            << std::endl;
-
   // post-smoothing
-  if (debug > 1)
-    deallog << cychar << "-cycle smoothing level " << level << std::endl;
   post_smooth->smooth(level, solution[level], defect2[level]);
-
-  if (debug > 2)
-    deallog << cychar << "-cycle solution norm   " << solution[level].l2_norm()
-            << std::endl;
-
-  if (debug > 1)
-    deallog << cychar << "-cycle leaving level   " << level << std::endl;
 }
 
 

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.