]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Deduce minimal and maximal valid Multigrid level from mg::Matrix 2948/head
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>
Thu, 15 Sep 2016 16:09:09 +0000 (18:09 +0200)
include/deal.II/multigrid/mg_base.h
include/deal.II/multigrid/mg_matrix.h
include/deal.II/multigrid/multigrid.h
include/deal.II/multigrid/multigrid.templates.h
tests/multigrid/cycles.cc
tests/multigrid/step-16-03.cc
tests/multigrid/step-16-04.cc

index 9c5ef476f4a5b36f498ba52164e2078d6de85e51..9129713cb642fa484b08b1ed715ed0cb4373459a 100644 (file)
@@ -79,6 +79,16 @@ public:
   virtual void Tvmult_add (const unsigned int level,
                            VectorType         &dst,
                            const VectorType   &src) const = 0;
+
+  /**
+   * Return the minimal level for which matrices are stored.
+   */
+  virtual unsigned int get_minlevel() const = 0;
+
+  /**
+   * Return the minimal level for which matrices are stored.
+   */
+  virtual unsigned int get_maxlevel() const = 0;
 };
 
 
index 66d1e8a33dc17d84eb9c93ea4275d0a31e5cffa7..0f6b257a63ab1ad20d679eac24504439b318fd96 100644 (file)
@@ -72,6 +72,8 @@ namespace mg
     virtual void vmult_add (const unsigned int level, VectorType &dst, const VectorType &src) const;
     virtual void Tvmult (const unsigned int level, VectorType &dst, const VectorType &src) const;
     virtual void Tvmult_add (const unsigned int level, VectorType &dst, const VectorType &src) const;
+    virtual unsigned int get_minlevel() const;
+    virtual unsigned int get_maxlevel() const;
 
     /**
      * Memory used by this object.
@@ -181,6 +183,7 @@ namespace mg
   }
 
 
+
   template <typename VectorType>
   template <typename MatrixType>
   inline
@@ -190,12 +193,14 @@ namespace mg
   }
 
 
+
   template <typename VectorType>
   inline
   Matrix<VectorType>::Matrix ()
   {}
 
 
+
   template <typename VectorType>
   inline
   const PointerMatrixBase<VectorType> &
@@ -205,6 +210,7 @@ namespace mg
   }
 
 
+
   template <typename VectorType>
   void
   Matrix<VectorType>::vmult (const unsigned int level,
@@ -215,6 +221,7 @@ namespace mg
   }
 
 
+
   template <typename VectorType>
   void
   Matrix<VectorType>::vmult_add (const unsigned int level,
@@ -225,6 +232,7 @@ namespace mg
   }
 
 
+
   template <typename VectorType>
   void
   Matrix<VectorType>::Tvmult (const unsigned int level,
@@ -235,6 +243,7 @@ namespace mg
   }
 
 
+
   template <typename VectorType>
   void
   Matrix<VectorType>::Tvmult_add (const unsigned int level,
@@ -245,6 +254,25 @@ namespace mg
   }
 
 
+
+  template <typename VectorType>
+  unsigned int
+  Matrix<VectorType>::get_minlevel() const
+  {
+    return matrices.min_level();
+  }
+
+
+
+  template <typename VectorType>
+  unsigned int
+  Matrix<VectorType>::get_maxlevel() const
+  {
+    return matrices.max_level();
+  }
+
+
+
   template <typename VectorType>
   inline
   std::size_t
@@ -278,6 +306,7 @@ MGMatrixSelect<MatrixType, number>::set_matrix (MGLevelObject<MatrixType> *p)
 }
 
 
+
 template <typename MatrixType, typename number>
 void
 MGMatrixSelect<MatrixType, number>::
@@ -289,6 +318,7 @@ select_block (const unsigned int brow,
 }
 
 
+
 template <typename MatrixType, typename number>
 void
 MGMatrixSelect<MatrixType, number>::
@@ -303,6 +333,7 @@ vmult  (const unsigned int    level,
 }
 
 
+
 template <typename MatrixType, typename number>
 void
 MGMatrixSelect<MatrixType, number>::
@@ -317,6 +348,7 @@ vmult_add  (const unsigned int    level,
 }
 
 
+
 template <typename MatrixType, typename number>
 void
 MGMatrixSelect<MatrixType, number>::
@@ -331,6 +363,7 @@ Tvmult  (const unsigned int    level,
 }
 
 
+
 template <typename MatrixType, typename number>
 void
 MGMatrixSelect<MatrixType, number>::
index 749b920a112ad2f5ba9015723ed05b70da460e99..f35cbf54747cf257a00a7d14edc287abec75841c 100644 (file)
@@ -91,6 +91,9 @@ public:
    * This function already initializes the vectors which will be used later in
    * the course of the computations. You should therefore create objects of
    * this type as late as possible.
+   *
+   * @deprecated Use the other constructor instead.
+   * The DoFHandler is actually not needed.
    */
   template <int dim>
   Multigrid(const DoFHandler<dim>              &mg_dof_handler,
@@ -99,20 +102,26 @@ public:
             const MGTransferBase<VectorType>   &transfer,
             const MGSmootherBase<VectorType>   &pre_smooth,
             const MGSmootherBase<VectorType>   &post_smooth,
-            Cycle                              cycle = v_cycle,
             const unsigned int                 minlevel = 0,
-            const unsigned int                 maxlevel = numbers::invalid_unsigned_int);
+            const unsigned int                 maxlevel = numbers::invalid_unsigned_int,
+            Cycle                              cycle = v_cycle) DEAL_II_DEPRECATED;
 
   /**
-   * Constructor. Same as above but without checking validity of minlevel and maxlevel.
+   * Constructor. <tt>transfer</tt> is an object performing prolongation and
+   * restriction. For levels in [minlevel, maxlevel] matrix has to contain
+   * valid matrices. By default the maxlevel is set to the maximal valid level.
+   *
+   * This function already initializes the vectors which will be used later in
+   * the course of the computations. You should therefore create objects of
+   * this type as late as possible.
    */
-  Multigrid(const unsigned int                 minlevel,
-            const unsigned int                 maxlevel,
-            const MGMatrixBase<VectorType>     &matrix,
+  Multigrid(const MGMatrixBase<VectorType>     &matrix,
             const MGCoarseGridBase<VectorType> &coarse,
             const MGTransferBase<VectorType>   &transfer,
             const MGSmootherBase<VectorType>   &pre_smooth,
             const MGSmootherBase<VectorType>   &post_smooth,
+            const unsigned int                 minlevel = 0,
+            const unsigned int                 maxlevel = numbers::invalid_unsigned_int,
             Cycle                              cycle = v_cycle);
 
   /**
@@ -432,9 +441,9 @@ Multigrid<VectorType>::Multigrid (const DoFHandler<dim>              &mg_dof_han
                                   const MGTransferBase<VectorType>   &transfer,
                                   const MGSmootherBase<VectorType>   &pre_smooth,
                                   const MGSmootherBase<VectorType>   &post_smooth,
-                                  Cycle                              cycle,
                                   const unsigned int                 min_level,
-                                  const unsigned int                 max_level)
+                                  const unsigned int                 max_level,
+                                  Cycle                              cycle)
   :
   cycle_type(cycle),
   minlevel(min_level),
@@ -454,15 +463,38 @@ Multigrid<VectorType>::Multigrid (const DoFHandler<dim>              &mg_dof_han
   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));
+  reinit(minlevel, maxlevel);
+}
+
+
 
-  defect.resize(minlevel,maxlevel);
-  solution.resize(minlevel,maxlevel);
-  t.resize(minlevel,maxlevel);
-  defect2.resize(minlevel,maxlevel);
+template <typename VectorType>
+Multigrid<VectorType>::Multigrid (const MGMatrixBase<VectorType>        &matrix,
+                                  const MGCoarseGridBase<VectorType>    &coarse,
+                                  const MGTransferBase<VectorType>      &transfer,
+                                  const MGSmootherBase<VectorType>      &pre_smooth,
+                                  const MGSmootherBase<VectorType>      &post_smooth,
+                                  const unsigned int                    min_level,
+                                  const unsigned int                    max_level,
+                                  Cycle                                 cycle)
+  :
+  cycle_type(cycle),
+  matrix(&matrix, typeid(*this).name()),
+  coarse(&coarse, typeid(*this).name()),
+  transfer(&transfer, typeid(*this).name()),
+  pre_smooth(&pre_smooth, typeid(*this).name()),
+  post_smooth(&post_smooth, typeid(*this).name()),
+  edge_out(0, typeid(*this).name()),
+  edge_in(0, typeid(*this).name()),
+  edge_down(0, typeid(*this).name()),
+  edge_up(0, typeid(*this).name()),
+  debug(0)
+{
+  if (max_level == numbers::invalid_unsigned_int)
+    maxlevel = matrix.get_maxlevel();
+  else
+    maxlevel = max_level;
+  reinit(min_level, maxlevel);
 }
 
 
index e185cd98cbc8ea88edeef8d68c579e3be045b2b9..7908b99a2f28580bdd1444c521b0e6d677afb5f2 100644 (file)
 
 DEAL_II_NAMESPACE_OPEN
 
-
-template <typename VectorType>
-Multigrid<VectorType>::Multigrid (const unsigned int                    minlevel,
-                                  const unsigned int                    maxlevel,
-                                  const MGMatrixBase<VectorType>        &matrix,
-                                  const MGCoarseGridBase<VectorType>    &coarse,
-                                  const MGTransferBase<VectorType>      &transfer,
-                                  const MGSmootherBase<VectorType>      &pre_smooth,
-                                  const MGSmootherBase<VectorType>      &post_smooth,
-                                  typename Multigrid<VectorType>::Cycle cycle)
-  :
-  cycle_type(cycle),
-  minlevel(minlevel),
-  maxlevel(maxlevel),
-  defect(minlevel,maxlevel),
-  solution(minlevel,maxlevel),
-  t(minlevel,maxlevel),
-  matrix(&matrix, typeid(*this).name()),
-  coarse(&coarse, typeid(*this).name()),
-  transfer(&transfer, typeid(*this).name()),
-  pre_smooth(&pre_smooth, typeid(*this).name()),
-  post_smooth(&post_smooth, typeid(*this).name()),
-  edge_out(0, typeid(*this).name()),
-  edge_in(0, typeid(*this).name()),
-  edge_down(0, typeid(*this).name()),
-  edge_up(0, typeid(*this).name()),
-  debug(0)
-{}
-
-
-
 template <typename VectorType>
 void
 Multigrid<VectorType>::reinit (const unsigned int min_level,
                                const unsigned int max_level)
 {
+  Assert (min_level >= matrix->get_minlevel(),
+          ExcLowerRangeType<unsigned int>(min_level, matrix->get_minlevel()));
+  Assert (max_level <= matrix->get_maxlevel(),
+          ExcLowerRangeType<unsigned int>(matrix->get_maxlevel(), max_level));
   Assert (min_level <= max_level,
           ExcLowerRangeType<unsigned int>(max_level, min_level));
   minlevel=min_level;
@@ -68,6 +41,7 @@ Multigrid<VectorType>::reinit (const unsigned int min_level,
 }
 
 
+
 template <typename VectorType>
 void
 Multigrid<VectorType>::set_maxlevel (const unsigned int l)
@@ -76,6 +50,7 @@ Multigrid<VectorType>::set_maxlevel (const unsigned int l)
 }
 
 
+
 template <typename VectorType>
 void
 Multigrid<VectorType>::set_minlevel (const unsigned int l,
@@ -88,6 +63,7 @@ Multigrid<VectorType>::set_minlevel (const unsigned int l,
 }
 
 
+
 template <typename VectorType>
 void
 Multigrid<VectorType>::set_cycle(typename Multigrid<VectorType>::Cycle c)
@@ -96,6 +72,7 @@ Multigrid<VectorType>::set_cycle(typename Multigrid<VectorType>::Cycle c)
 }
 
 
+
 template <typename VectorType>
 void
 Multigrid<VectorType>::set_debug (const unsigned int d)
@@ -104,6 +81,7 @@ Multigrid<VectorType>::set_debug (const unsigned int d)
 }
 
 
+
 template <typename VectorType>
 void
 Multigrid<VectorType>::set_edge_matrices (const MGMatrixBase<VectorType> &down,
@@ -114,6 +92,7 @@ Multigrid<VectorType>::set_edge_matrices (const MGMatrixBase<VectorType> &down,
 }
 
 
+
 template <typename VectorType>
 void
 Multigrid<VectorType>::set_edge_flux_matrices (const MGMatrixBase<VectorType> &down,
@@ -124,6 +103,7 @@ Multigrid<VectorType>::set_edge_flux_matrices (const MGMatrixBase<VectorType> &d
 }
 
 
+
 template <typename VectorType>
 void
 Multigrid<VectorType>::level_v_step (const unsigned int level)
@@ -374,6 +354,7 @@ Multigrid<VectorType>::level_step(const unsigned int level,
 }
 
 
+
 template <typename VectorType>
 void
 Multigrid<VectorType>::cycle()
@@ -401,6 +382,7 @@ Multigrid<VectorType>::cycle()
 }
 
 
+
 template <typename VectorType>
 void
 Multigrid<VectorType>::vcycle()
index 05e607f9d069a2503f732e5fd967d8d4907e7681..53edc311ecb16d8ce968aecdd304467c4cf24f41 100644 (file)
@@ -72,8 +72,8 @@ void test_cycles(unsigned int minlevel, unsigned int maxlevel)
     level_matrices[i].reinit(N, N);
   mg::Matrix<VectorType> mgmatrix(level_matrices);
 
-  Multigrid<VectorType> mg1(minlevel, maxlevel, mgmatrix, all, all, all, all,
-                            Multigrid<VectorType>::v_cycle);
+  Multigrid<VectorType> mg1(mgmatrix, all, all, all, all,
+                            minlevel, maxlevel, Multigrid<VectorType>::v_cycle);
   mg1.set_debug(3);
   for (unsigned int i=minlevel; i<=maxlevel; ++i)
     mg1.defect[i].reinit(N);
index 0d069340bd2f19e09d8861df77139c4ba7e3b339..63ac2e4b1a3051fcb61c561bd354fdd5a24e337f 100644 (file)
@@ -390,9 +390,9 @@ void LaplaceProblem<dim>::solve ()
                                 mg_transfer,
                                 mg_smoother,
                                 mg_smoother,
-                                Multigrid<Vector<double> >::v_cycle,
                                 min_level,
-                                triangulation.n_global_levels()-1);
+                                triangulation.n_global_levels()-1,
+                                Multigrid<Vector<double> >::v_cycle);
   mg.set_edge_matrices(mg_interface_down, mg_interface_up);
 
   PreconditionMG<dim, Vector<double>, MGTransferPrebuilt<Vector<double> > >
index 65b4d650716fd3edfe975377ec21bd8cd14e579d..4e70d17b58f6ea228b4641c1550b1abecb86c347 100644 (file)
@@ -396,13 +396,15 @@ void LaplaceProblem<dim>::solve ()
   mg::Matrix<LinearAlgebra::distributed::Vector<double> > mg_interface_up(mg_interface_matrices);
   mg::Matrix<LinearAlgebra::distributed::Vector<double> > mg_interface_down(mg_interface_matrices);
 
-  Multigrid<LinearAlgebra::distributed::Vector<double> > mg(min_level,
-                                                            triangulation.n_global_levels()-1,
-                                                            mg_matrix,
+  Multigrid<LinearAlgebra::distributed::Vector<double> > mg(mg_matrix,
                                                             coarse_grid_solver,
                                                             mg_transfer,
                                                             mg_smoother,
-                                                            mg_smoother);
+                                                            mg_smoother,
+                                                            min_level,
+                                                            triangulation.n_global_levels()-1);
+  Assert(min_level == mg.get_minlevel(), ExcInternalError());
+  Assert(triangulation.n_global_levels()-1 == mg.get_maxlevel(), ExcInternalError());
   mg.set_edge_matrices(mg_interface_down, mg_interface_up);
 
   PreconditionMG<dim, LinearAlgebra::distributed::Vector<double>, MGTransferPrebuilt<LinearAlgebra::distributed::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.