From 03ebce802508802d6d421c028c20b70e6cc88ebb Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Thu, 11 Jul 2019 16:56:22 -0400 Subject: [PATCH] Remove deprecated variables in PreconditionChebyshev::AdditionalData --- include/deal.II/lac/precondition.h | 49 ++-------------------- tests/lac/precondition_chebyshev_01.cc | 6 ++- tests/matrix_free/multigrid_dg_periodic.cc | 6 ++- tests/matrix_free/multigrid_dg_sip_01.cc | 6 ++- tests/matrix_free/multigrid_dg_sip_02.cc | 6 ++- tests/matrix_free/parallel_multigrid_mf.cc | 6 ++- 6 files changed, 24 insertions(+), 55 deletions(-) diff --git a/include/deal.II/lac/precondition.h b/include/deal.II/lac/precondition.h index b840a8087d..52a518b141 100644 --- a/include/deal.II/lac/precondition.h +++ b/include/deal.II/lac/precondition.h @@ -1024,7 +1024,6 @@ public: */ AdditionalData(const unsigned int degree = 1, const double smoothing_range = 0., - const bool nonzero_starting = false, const unsigned int eig_cg_n_iterations = 8, const double eig_cg_residual = 1e-2, const double max_eigenvalue = 1); @@ -1055,20 +1054,6 @@ public: */ double smoothing_range; - /** - * When this flag is set to true, it enables the method - * vmult(dst, src) to use non-zero data in the vector - * dst, appending to it the Chebyshev corrections. This can be - * useful in some situations (e.g. when used for high-frequency error - * smoothing in a multigrid algorithm), but not the way the solver classes - * expect a preconditioner to work (where one ignores the content in - * dst for the preconditioner application). - * - * @deprecated For non-zero starting, use the step() and Tstep() - * interfaces, whereas vmult() provides the preconditioner interface. - */ - bool nonzero_starting DEAL_II_DEPRECATED; - /** * Maximum number of CG iterations performed for finding the maximum * eigenvalue. If set to zero, no computations are performed. Instead, the @@ -1090,13 +1075,6 @@ public: */ double max_eigenvalue; - /** - * Stores the inverse of the diagonal of the underlying matrix. - * - * @deprecated Set the variable @p preconditioner defined below instead. - */ - VectorType matrix_diagonal_inverse DEAL_II_DEPRECATED; - /** * Stores the preconditioner object that the Chebyshev is wrapped around. */ @@ -2008,14 +1986,11 @@ namespace internal solution.swap(solution_old); } - template + template inline void initialize_preconditioner( const MatrixType & matrix, - std::shared_ptr &preconditioner, - VectorType &) + std::shared_ptr &preconditioner) { (void)matrix; (void)preconditioner; @@ -2026,8 +2001,7 @@ namespace internal inline void initialize_preconditioner( const MatrixType & matrix, - std::shared_ptr> &preconditioner, - VectorType & diagonal_inverse) + std::shared_ptr> &preconditioner) { if (preconditioner.get() == nullptr || preconditioner->m() != matrix.m()) { @@ -2039,14 +2013,6 @@ namespace internal ExcMessage( "Preconditioner appears to be initialized but not sized correctly")); - // Check if we can initialize from vector that then gets set to zero - // as the matrix will own the memory - preconditioner->reinit(diagonal_inverse); - { - VectorType empty_vector; - diagonal_inverse.reinit(empty_vector); - } - // This part only works in serial if (preconditioner->m() != matrix.m()) { @@ -2118,19 +2084,15 @@ namespace internal -// avoid warning about deprecated variable nonzero_starting - template inline PreconditionChebyshev:: AdditionalData::AdditionalData(const unsigned int degree, const double smoothing_range, - const bool nonzero_starting, const unsigned int eig_cg_n_iterations, const double eig_cg_residual, const double max_eigenvalue) : degree(degree) , smoothing_range(smoothing_range) - , nonzero_starting(nonzero_starting) , eig_cg_n_iterations(eig_cg_n_iterations) , eig_cg_residual(eig_cg_residual) , max_eigenvalue(max_eigenvalue) @@ -2151,8 +2113,6 @@ inline PreconditionChebyshev:: } -// avoid warning about deprecated variable -// AdditionalData::matrix_diagonal_inverse template inline void @@ -2165,7 +2125,7 @@ PreconditionChebyshev::initialize( Assert(data.degree > 0, ExcMessage("The degree of the Chebyshev method must be positive.")); internal::PreconditionChebyshevImplementation::initialize_preconditioner( - matrix, data.preconditioner, data.matrix_diagonal_inverse); + matrix, data.preconditioner); eigenvalues_are_initialized = false; } @@ -2180,7 +2140,6 @@ PreconditionChebyshev::clear() matrix_ptr = nullptr; { VectorType empty_vector; - data.matrix_diagonal_inverse.reinit(empty_vector); solution_old.reinit(empty_vector); temp_vector1.reinit(empty_vector); temp_vector2.reinit(empty_vector); diff --git a/tests/lac/precondition_chebyshev_01.cc b/tests/lac/precondition_chebyshev_01.cc index 929065b16b..7af71d4388 100644 --- a/tests/lac/precondition_chebyshev_01.cc +++ b/tests/lac/precondition_chebyshev_01.cc @@ -77,8 +77,10 @@ check() deallog << std::endl; Vector matrix_diagonal(size); - matrix_diagonal = 1; - data.matrix_diagonal_inverse = matrix_diagonal; + matrix_diagonal = 1; + auto preconditioner = std::make_shared>>(); + preconditioner->reinit(matrix_diagonal); + data.preconditioner = std::move(preconditioner); prec.initialize(m, data); deallog << "Check vmult diag: "; diff --git a/tests/matrix_free/multigrid_dg_periodic.cc b/tests/matrix_free/multigrid_dg_periodic.cc index 79199c7c12..853d5e0fe4 100644 --- a/tests/matrix_free/multigrid_dg_periodic.cc +++ b/tests/matrix_free/multigrid_dg_periodic.cc @@ -572,8 +572,10 @@ do_test(const DoFHandler &dof) smoother_data[level].smoothing_range = 20.; smoother_data[level].degree = 5; smoother_data[level].eig_cg_n_iterations = 15; - smoother_data[level].matrix_diagonal_inverse = - mg_matrices[level].get_matrix_diagonal_inverse(); + auto preconditioner = std::make_shared< + DiagonalMatrix>>(); + preconditioner->reinit(mg_matrices[level].get_matrix_diagonal_inverse()); + smoother_data[level].preconditioner = std::move(preconditioner); } mg_smoother.initialize(mg_matrices, smoother_data); diff --git a/tests/matrix_free/multigrid_dg_sip_01.cc b/tests/matrix_free/multigrid_dg_sip_01.cc index 5c65da9159..23096ad0bf 100644 --- a/tests/matrix_free/multigrid_dg_sip_01.cc +++ b/tests/matrix_free/multigrid_dg_sip_01.cc @@ -582,8 +582,10 @@ do_test(const DoFHandler &dof, const bool also_test_parallel = false) smoother_data[level].smoothing_range = 20.; smoother_data[level].degree = 5; smoother_data[level].eig_cg_n_iterations = 15; - smoother_data[level].matrix_diagonal_inverse = - mg_matrices[level].get_matrix_diagonal_inverse(); + auto preconditioner = std::make_shared< + DiagonalMatrix>>(); + preconditioner->reinit(mg_matrices[level].get_matrix_diagonal_inverse()); + smoother_data[level].preconditioner = std::move(preconditioner); } mg_smoother.initialize(mg_matrices, smoother_data); diff --git a/tests/matrix_free/multigrid_dg_sip_02.cc b/tests/matrix_free/multigrid_dg_sip_02.cc index 3afcea7283..52181a9334 100644 --- a/tests/matrix_free/multigrid_dg_sip_02.cc +++ b/tests/matrix_free/multigrid_dg_sip_02.cc @@ -494,8 +494,10 @@ do_test(const DoFHandler &dof, const bool also_test_parallel = false) smoother_data[level].smoothing_range = 20.; smoother_data[level].degree = 5; smoother_data[level].eig_cg_n_iterations = 15; - smoother_data[level].matrix_diagonal_inverse = - mg_matrices[level].get_matrix_diagonal_inverse(); + auto preconditioner = std::make_shared< + DiagonalMatrix>>(); + preconditioner->reinit(mg_matrices[level].get_matrix_diagonal_inverse()); + smoother_data[level].preconditioner = std::move(preconditioner); } mg_smoother.initialize(mg_matrices, smoother_data); diff --git a/tests/matrix_free/parallel_multigrid_mf.cc b/tests/matrix_free/parallel_multigrid_mf.cc index 44952edee9..5d4dc67743 100644 --- a/tests/matrix_free/parallel_multigrid_mf.cc +++ b/tests/matrix_free/parallel_multigrid_mf.cc @@ -405,8 +405,10 @@ do_test(const DoFHandler &dof) smoother_data[level].smoothing_range = 15.; smoother_data[level].degree = 5; smoother_data[level].eig_cg_n_iterations = 15; - smoother_data[level].matrix_diagonal_inverse = - mg_matrices[level].get_matrix_diagonal_inverse(); + auto preconditioner = std::make_shared< + DiagonalMatrix>>(); + preconditioner->reinit(mg_matrices[level].get_matrix_diagonal_inverse()); + smoother_data[level].preconditioner = std::move(preconditioner); } // temporarily disable deallog for the setup of the preconditioner that -- 2.39.5