From: Daniel Arndt <daniel.arndt@iwr.uni-heidelberg.de> Date: Mon, 15 Aug 2016 12:04:38 +0000 (+0200) Subject: Deduce minimal and maximal valid Multigrid level from mg::Matrix X-Git-Tag: v8.5.0-rc1~642^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=3947a1c22e37bcc50755c6705986246aedcb2004;p=dealii.git Deduce minimal and maximal valid Multigrid level from mg::Matrix --- diff --git a/include/deal.II/multigrid/mg_base.h b/include/deal.II/multigrid/mg_base.h index 9c5ef476f4..9129713cb6 100644 --- a/include/deal.II/multigrid/mg_base.h +++ b/include/deal.II/multigrid/mg_base.h @@ -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; }; diff --git a/include/deal.II/multigrid/mg_matrix.h b/include/deal.II/multigrid/mg_matrix.h index 66d1e8a33d..0f6b257a63 100644 --- a/include/deal.II/multigrid/mg_matrix.h +++ b/include/deal.II/multigrid/mg_matrix.h @@ -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>:: diff --git a/include/deal.II/multigrid/multigrid.h b/include/deal.II/multigrid/multigrid.h index 749b920a11..f35cbf5474 100644 --- a/include/deal.II/multigrid/multigrid.h +++ b/include/deal.II/multigrid/multigrid.h @@ -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); } diff --git a/include/deal.II/multigrid/multigrid.templates.h b/include/deal.II/multigrid/multigrid.templates.h index e185cd98cb..7908b99a2f 100644 --- a/include/deal.II/multigrid/multigrid.templates.h +++ b/include/deal.II/multigrid/multigrid.templates.h @@ -23,42 +23,15 @@ 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() diff --git a/tests/multigrid/cycles.cc b/tests/multigrid/cycles.cc index 05e607f9d0..53edc311ec 100644 --- a/tests/multigrid/cycles.cc +++ b/tests/multigrid/cycles.cc @@ -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); diff --git a/tests/multigrid/step-16-03.cc b/tests/multigrid/step-16-03.cc index 0d069340bd..63ac2e4b1a 100644 --- a/tests/multigrid/step-16-03.cc +++ b/tests/multigrid/step-16-03.cc @@ -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> > > diff --git a/tests/multigrid/step-16-04.cc b/tests/multigrid/step-16-04.cc index 65b4d65071..4e70d17b58 100644 --- a/tests/multigrid/step-16-04.cc +++ b/tests/multigrid/step-16-04.cc @@ -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> > >