]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Allow to specify minlevel and maxlevel in the Multigrid constructor
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Mon, 15 Aug 2016 12:04:38 +0000 (14:04 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Mon, 15 Aug 2016 12:33:20 +0000 (14:33 +0200)
include/deal.II/multigrid/multigrid.h
include/deal.II/multigrid/multigrid.templates.h
tests/multigrid/step-16-03.cc

index b8b4fa788a60438ec0179b8bba18d3b3ffff0bc1..749b920a112ad2f5ba9015723ed05b70da460e99 100644 (file)
@@ -81,8 +81,11 @@ public:
   typedef const VectorType const_vector_type;
 
   /**
-   * Constructor. The DoFHandler is used to determine the highest possible
-   * level. <tt>transfer</tt> is an object performing prolongation and
+   * Constructor. The DoFHandler is used to check whether the provided
+   * minlevel and maxlevel are in the range of valid levels.
+   * If maxlevel is set to the default value, the highest valid
+   * level is used.
+   * <tt>transfer</tt> is an object performing prolongation and
    * restriction.
    *
    * This function already initializes the vectors which will be used later in
@@ -96,11 +99,12 @@ public:
             const MGTransferBase<VectorType>   &transfer,
             const MGSmootherBase<VectorType>   &pre_smooth,
             const MGSmootherBase<VectorType>   &post_smooth,
-            Cycle                              cycle = v_cycle);
+            Cycle                              cycle = v_cycle,
+            const unsigned int                 minlevel = 0,
+            const unsigned int                 maxlevel = numbers::invalid_unsigned_int);
 
   /**
-   * Constructor. Same as above but you can determine the multigrid
-   * levels to use yourself.
+   * Constructor. Same as above but without checking validity of minlevel and maxlevel.
    */
   Multigrid(const unsigned int                 minlevel,
             const unsigned int                 maxlevel,
@@ -422,21 +426,18 @@ private:
 
 template <typename VectorType>
 template <int dim>
-Multigrid<VectorType>::Multigrid (const DoFHandler<dim>          &mg_dof_handler,
+Multigrid<VectorType>::Multigrid (const DoFHandler<dim>              &mg_dof_handler,
                                   const MGMatrixBase<VectorType>     &matrix,
                                   const MGCoarseGridBase<VectorType> &coarse,
                                   const MGTransferBase<VectorType>   &transfer,
                                   const MGSmootherBase<VectorType>   &pre_smooth,
                                   const MGSmootherBase<VectorType>   &post_smooth,
-                                  Cycle                              cycle)
+                                  Cycle                              cycle,
+                                  const unsigned int                 min_level,
+                                  const unsigned int                 max_level)
   :
   cycle_type(cycle),
-  minlevel(0),
-  maxlevel(mg_dof_handler.get_triangulation().n_global_levels()-1),
-  defect(minlevel,maxlevel),
-  solution(minlevel,maxlevel),
-  t(minlevel,maxlevel),
-  defect2(minlevel,maxlevel),
+  minlevel(min_level),
   matrix(&matrix, typeid(*this).name()),
   coarse(&coarse, typeid(*this).name()),
   transfer(&transfer, typeid(*this).name()),
@@ -445,7 +446,24 @@ Multigrid<VectorType>::Multigrid (const DoFHandler<dim>          &mg_dof_handler
   edge_down(0, typeid(*this).name()),
   edge_up(0, typeid(*this).name()),
   debug(0)
-{}
+{
+  const unsigned int dof_handler_max_level
+    = mg_dof_handler.get_triangulation().n_global_levels()-1;
+  if (max_level == numbers::invalid_unsigned_int)
+    maxlevel = dof_handler_max_level;
+  else
+    maxlevel = max_level;
+
+  Assert (maxlevel <= dof_handler_max_level,
+          ExcLowerRangeType<unsigned int>(dof_handler_max_level, max_level));
+  Assert (minlevel <= maxlevel,
+          ExcLowerRangeType<unsigned int>(max_level, min_level));
+
+  defect.resize(minlevel,maxlevel);
+  solution.resize(minlevel,maxlevel);
+  t.resize(minlevel,maxlevel);
+  defect2.resize(minlevel,maxlevel);
+}
 
 
 
index 5270febce77694a714c5bbb1e5b54425a4baf21d..e185cd98cbc8ea88edeef8d68c579e3be045b2b9 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 1999 - 2015 by the deal.II authors
+// Copyright (C) 1999 - 2016 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
index 4630093df1ffd15f4b92f71972d3e1d6f5ba0e54..0d069340bd2f19e09d8861df77139c4ba7e3b339 100644 (file)
@@ -384,13 +384,15 @@ void LaplaceProblem<dim>::solve ()
   mg::Matrix<> mg_interface_up(mg_interface_matrices);
   mg::Matrix<> mg_interface_down(mg_interface_matrices);
 
-  Multigrid<Vector<double> > mg(min_level,
-                                triangulation.n_global_levels()-1,
+  Multigrid<Vector<double> > mg(mg_dof_handler,
                                 mg_matrix,
                                 coarse_grid_solver,
                                 mg_transfer,
                                 mg_smoother,
-                                mg_smoother);
+                                mg_smoother,
+                                Multigrid<Vector<double> >::v_cycle,
+                                min_level,
+                                triangulation.n_global_levels()-1);
   mg.set_edge_matrices(mg_interface_down, mg_interface_up);
 
   PreconditionMG<dim, Vector<double>, MGTransferPrebuilt<Vector<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.