]> https://gitweb.dealii.org/ - dealii.git/commitdiff
W- and F-cycles running
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Sat, 19 Mar 2005 05:09:29 +0000 (05:09 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Sat, 19 Mar 2005 05:09:29 +0000 (05:09 +0000)
git-svn-id: https://svn.dealii.org/trunk@10190 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/multigrid/multigrid.h
deal.II/deal.II/include/multigrid/multigrid.templates.h

index 4b5ef4d7079472e1f69f110ba32fb39170f43f33..794f8534ea9f39b67ce86175465d58e2d548b91c 100644 (file)
@@ -40,7 +40,8 @@
  * The function which starts a multigrid cycle on the finest level is
  * cycle(). Depending on the cycle type chosen with the constructor
  * (see enum Cycle), this function triggers one of the cycles
- * level_v_step(), level_w_step() or level_f_step().
+ * level_v_step() or level_step(), where the latter one can do
+ * different types of cycles.
  *
  * @note The interface of this class is still very clumsy. In
  * particular, you will have to set up quite a few auxiliary objects
@@ -247,18 +248,20 @@ class Multigrid : public Subscriptor
     void level_v_step (const unsigned int level);
 
                                     /**
-                                     * The actual W-cycle multigrid method.
-                                     * <tt>level</tt> is the level the
-                                     * function starts on. It
+                                     * The actual W-cycle or F-cycle
+                                     * multigrid method.
+                                     * <tt>level</tt> is the level
+                                     * the function starts on. It
                                      * will usually be called for the
                                      * highest level from outside,
                                      * but will then call itself
-                                     * recursively for <tt>level-1</tt>,
-                                     * unless we are on #minlevel
-                                     * where the coarse grid solver
-                                     * solves the problem exactly.
+                                     * recursively for
+                                     * <tt>level-1</tt>, unless we
+                                     * are on #minlevel where the
+                                     * coarse grid solver solves the
+                                     * problem exactly.
                                      */
-    void level_w_step (const unsigned int level);
+    void level_step (const unsigned int level, Cycle cycle);
 
                                     /**
                                      * Cycle type performed by the method cycle().
index 2c36f0114e4be7c28158e0c297cfd4e7f8472370..5cacf07327676dfcf437baf8b925298b16b3f278 100644 (file)
 using namespace std;
 
 
+
+
+
+
+
+
+
 template <class VECTOR>
 Multigrid<VECTOR>::Multigrid (const unsigned int minlevel,
                              const unsigned int maxlevel,
@@ -188,14 +195,31 @@ Multigrid<VECTOR>::level_v_step(const unsigned int level)
 
 template <class VECTOR>
 void
-Multigrid<VECTOR>::level_w_step(const unsigned int level)
+Multigrid<VECTOR>::level_step(const unsigned int level,
+                             Cycle cycle)
 {
+                                  // Not actually the defect yet, but
+                                  // we do not want to spend yet
+                                  // another vector.
+  if (level>minlevel)
+    {
+      defect2[level-1] = 0.;
+      transfer->restrict_and_add (level, defect2[level-1], defect2[level]);
+    }
+  
+                                  // We get an update of the defect
+                                  // from the previous level in t and
+                                  // from two levels above in
+                                  // defect2. This is subtracted from
+                                  // the original defect.
+  t[level].equ(-1.,defect2[level],1.,defect[level]);
+  
   if (level == minlevel)
     {
       if (debug>0)
        deallog << "W-cycle solving  level " << level << std::endl;
   
-      (*coarse)(level, solution[level], defect[level]);
+      (*coarse)(level, solution[level], t[level]);
       return;
     }
   if (debug>0)
@@ -205,55 +229,57 @@ Multigrid<VECTOR>::level_w_step(const unsigned int level)
     deallog << "Smoothing on     level " << level << std::endl;
                                   // smoothing of the residual by
                                   // modifying s
-  pre_smooth->smooth(level, solution[level], defect[level]);
+  pre_smooth->smooth(level, solution[level], t[level]);
 
   if (debug>1)
     deallog << "Residual on      level " << level << std::endl;
                                   // 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.
-  for (unsigned int l = level;l>minlevel;--l)
-    {
-      t[l-1] = 0.;
-      if (l==level && edge_down != 0)
-       {
-         edge_down->vmult(level, t[level-1], solution[level]);
-       }
-      transfer->restrict_and_add (l, t[l-1], t[l]);
-      defect[l-1] -= t[l-1];
-    }
-
+  if (edge_down != 0)
+    edge_down->vmult_add(level, defect2[level-1], solution[level]);
+  transfer->restrict_and_add (level, defect2[level-1], t[level]);
                                   // do recursion
   solution[level-1] = 0.;
-  level_w_step(level-1);
-  level_w_step(level-1);
+                                  // Every cycle need one recursion,
+                                  // the V-cycle, which is included
+                                  // here for the sake of the
+                                  // F-cycle, needs only one,
+  level_step(level-1, cycle);
+                                  // while the W-cycle repeats itself
+  if (cycle == w_cycle)
+    level_step(level-1, cycle);
+                                  // and the F-cycle does a V-cycle
+                                  // after an F-cycle.
+  else if (cycle == f_cycle)
+    level_step(level-1, v_cycle);
 
                                   // reset size of the auxiliary
                                   // vector, since it has been
                                   // resized in the recursive call to
                                   // level_v_step directly above
   t[level] = 0.;
-
                                   // do coarse grid correction
   transfer->prolongate(level, t[level], solution[level-1]);
   
   solution[level] += t[level];
 
+  
   if (edge_up != 0)
     {
       edge_up->Tvmult(level, t[level], solution[level-1]);
-      defect[level] -= t[level];
     }
-  
+
+  t[level].sadd(-1.,-1.,defect2[level],1.,defect[level]);
+
   if (debug>1)
     deallog << "Smoothing on     level " << level << std::endl;
                                   // post-smoothing
-  post_smooth->smooth(level, solution[level], defect[level]);
+  post_smooth->smooth(level, solution[level], t[level]);
 
   if (debug>1)
     deallog << "W-cycle leaving  level " << level << std::endl;
@@ -269,24 +295,21 @@ Multigrid<VECTOR>::cycle()
                                   // adjust the other vectors.
   solution.resize(minlevel, maxlevel);
   t.resize(minlevel, maxlevel);
-  
+  if (cycle_type != v_cycle)
+    defect2.resize(minlevel, maxlevel);
+    
   for (unsigned int level=minlevel; level<=maxlevel;++level)
     {
       solution[level].reinit(defect[level]);
       t[level].reinit(defect[level]);
+      if (cycle_type != v_cycle)
+       defect2[level].reinit(defect[level]);
     }
   
-  switch(cycle_type)
-    {
-      case v_cycle:
-       level_v_step (maxlevel);
-       break;
-      case w_cycle:
-       level_w_step (maxlevel);
-       break;
-      case f_cycle:
-       break;
-    }
+  if (cycle_type == v_cycle)
+    level_v_step (maxlevel);
+  else
+    level_step (maxlevel, cycle_type);
 }
 
 

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.