From 11b60c16e5f2cc0c3eb3ab9be2d8b96c9bb500d6 Mon Sep 17 00:00:00 2001 From: David Wells Date: Wed, 27 May 2020 15:39:45 -0400 Subject: [PATCH] Deprecate some old options for SparseVanka. --- include/deal.II/lac/sparse_vanka.h | 77 ++++++++--------- include/deal.II/lac/sparse_vanka.templates.h | 91 +++++++++----------- 2 files changed, 76 insertions(+), 92 deletions(-) diff --git a/include/deal.II/lac/sparse_vanka.h b/include/deal.II/lac/sparse_vanka.h index 07c5c2db91..dbd9c81ecf 100644 --- a/include/deal.II/lac/sparse_vanka.h +++ b/include/deal.II/lac/sparse_vanka.h @@ -153,6 +153,18 @@ public: */ SparseVanka(); + /** + * Constructor which also takes two deprecated inputs. + * + * @deprecated The use of the last two parameters is deprecated. They are + * currently ignored. + */ + DEAL_II_DEPRECATED + SparseVanka(const SparseMatrix &M, + const std::vector & selected, + const bool conserve_memory, + const unsigned int n_threads = MultithreadInfo::n_threads()); + /** * Constructor. Gets the matrix for preconditioning and a bit vector with * entries @p true for all rows to be updated. A reference to this vector @@ -164,24 +176,8 @@ public: * conceivable that the preconditioner is build up for one matrix once, but * is used for subsequent steps in a nonlinear process as well, where the * matrix changes in each step slightly. - * - * If @p conserve_mem is @p false, then the inverses of the local systems - * are computed now; if the flag is @p true, then they are computed every - * time the preconditioner is applied. This saves some memory, but makes - * preconditioning very slow. Note also, that if the flag is @p false, then - * the contents of the matrix @p M at the time of calling this constructor - * are used, while if the flag is @p true, then the values in @p M at the - * time of preconditioning are used. This may lead to different results, - * obviously, of @p M changes. - * - * The parameter @p n_threads determines how many threads shall be used in - * parallel when building the inverses of the diagonal blocks. This - * parameter is ignored if not in multithreaded mode. */ - SparseVanka(const SparseMatrix &M, - const std::vector & selected, - const bool conserve_memory = false, - const unsigned int n_threads = MultithreadInfo::n_threads()); + SparseVanka(const SparseMatrix &M, const std::vector &selected); /** * Destructor. Delete all allocated matrices. @@ -197,25 +193,23 @@ public: /** * Constructor. For the parameters' description, see below. */ + explicit AdditionalData(const std::vector &selected); + + /** + * Constructor. For the parameters' description, see below. + * + * @deprecated The use of this constructor is deprecated - the second and + * third parameters are ignored. + */ + DEAL_II_DEPRECATED AdditionalData(const std::vector &selected, - const bool conserve_memory = false, + const bool conserve_memory, const unsigned int n_threads = MultithreadInfo::n_threads()); /** * Indices of those degrees of freedom that we shall work on. */ const std::vector &selected; - - /** - * Conserve memory flag. - */ - const bool conserve_mem; - - /** - * Number of threads to be used when building the inverses. Only relevant - * in multithreaded mode. - */ - const unsigned int n_threads; }; @@ -310,22 +304,11 @@ private: */ SmartPointer, SparseVanka> matrix; - /** - * Conserve memory flag. - */ - bool conserve_mem; - /** * Indices of those degrees of freedom that we shall work on. */ const std::vector *selected; - /** - * Number of threads to be used when building the inverses. Only relevant in - * multithreaded mode. - */ - unsigned int n_threads; - /** * Array of inverse matrices, one for each degree of freedom. Only those * elements will be used that are tagged in @p selected. @@ -544,14 +527,26 @@ public: /** * Constructor. Pass all arguments except for @p n_blocks to the base class. + * + * @deprecated This constructor is deprecated. The values passed to the last + * two arguments are ignored. */ + DEAL_II_DEPRECATED SparseBlockVanka(const SparseMatrix &M, const std::vector & selected, const unsigned int n_blocks, const BlockingStrategy blocking_strategy, - const bool conserve_memory = false, + const bool conserve_memory, const unsigned int n_threads = MultithreadInfo::n_threads()); + /** + * Constructor. Pass all arguments except for @p n_blocks to the base class. + */ + SparseBlockVanka(const SparseMatrix &M, + const std::vector & selected, + const unsigned int n_blocks, + const BlockingStrategy blocking_strategy); + /** * Apply the preconditioner. */ diff --git a/include/deal.II/lac/sparse_vanka.templates.h b/include/deal.II/lac/sparse_vanka.templates.h index 823c670c2c..00fc4937d7 100644 --- a/include/deal.II/lac/sparse_vanka.templates.h +++ b/include/deal.II/lac/sparse_vanka.templates.h @@ -35,9 +35,7 @@ DEAL_II_NAMESPACE_OPEN template SparseVanka::SparseVanka() : matrix() - , conserve_mem(false) , selected() - , n_threads(0) , inverses() , _m(0) , _n(0) @@ -46,12 +44,16 @@ SparseVanka::SparseVanka() template SparseVanka::SparseVanka(const SparseMatrix &M, const std::vector & selected_dofs, - const bool conserve_mem, - const unsigned int n_threads) + const bool /*conserve_mem*/, + const unsigned int /*n_threads*/) + : SparseVanka(M, selected_dofs) +{} + +template +SparseVanka::SparseVanka(const SparseMatrix &M, + const std::vector & selected_dofs) : matrix(&M, typeid(*this).name()) - , conserve_mem(conserve_mem) , selected(&selected_dofs) - , n_threads(n_threads) , inverses(M.m(), nullptr) , _m(M.m()) , _n(M.n()) @@ -60,8 +62,7 @@ SparseVanka::SparseVanka(const SparseMatrix &M, Assert(M.m() == selected->size(), ExcDimensionMismatch(M.m(), selected->size())); - if (conserve_mem == false) - compute_inverses(); + compute_inverses(); } @@ -85,10 +86,8 @@ void SparseVanka::initialize(const SparseMatrix &M, const AdditionalData & additional_data) { - matrix = &M; - conserve_mem = additional_data.conserve_mem; - selected = &(additional_data.selected); - n_threads = additional_data.n_threads; + matrix = &M; + selected = &(additional_data.selected); inverses.resize(M.m()); _m = M.m(); _n = M.n(); @@ -97,8 +96,7 @@ SparseVanka::initialize(const SparseMatrix &M, Assert(M.m() == selected->size(), ExcDimensionMismatch(M.m(), selected->size())); - if (conserve_mem == false) - compute_inverses(); + compute_inverses(); } template @@ -113,8 +111,8 @@ SparseVanka::compute_inverses() #else const size_type n_inverses = std::count(selected->begin(), selected->end(), true); - - const size_type n_inverses_per_thread = + const std::size_t n_threads = MultithreadInfo::n_threads(); + const size_type n_inverses_per_thread = std::max(n_inverses / n_threads, static_cast(1U)); // set up start and end index @@ -285,17 +283,6 @@ SparseVanka::apply_preconditioner( { const size_type row_length = structure.row_length(row); - // if we don't store the - // inverse matrices, then alias - // the entry in the global - // vector to the local matrix - // to be used - if (conserve_mem == true) - { - inverses[row] = &local_matrix; - inverses[row]->reinit(row_length, row_length); - } - b.reinit(row_length); x.reinit(row_length); // mapping between: @@ -361,18 +348,9 @@ SparseVanka::apply_preconditioner( ((*dof_mask)[p->column()] == true)) b(i) -= p->value() * dst(p->column()); } - else - // if so, then build the - // matrix out of it - if (conserve_mem == true) - (*inverses[row])(i, js->second) = p->value(); } } - // Compute new values - if (conserve_mem == true) - inverses[row]->gauss_jordan(); - // apply preconditioner inverses[row]->vmult(x, b); @@ -390,12 +368,6 @@ SparseVanka::apply_preconditioner( // do nothing if not in // the range } - - // if we don't store the - // inverses, then unalias the - // local matrix - if (conserve_mem == true) - inverses[row] = nullptr; } } @@ -414,14 +386,21 @@ SparseVanka::memory_consumption() const } + template SparseVanka::AdditionalData::AdditionalData( - const std::vector &selected, - const bool conserve_mem, - const unsigned int n_threads) + const std::vector &selected) : selected(selected) - , conserve_mem(conserve_mem) - , n_threads(n_threads) +{} + + + +template +SparseVanka::AdditionalData::AdditionalData( + const std::vector &selected, + const bool /*conserve_mem*/, + const unsigned int /*n_threads*/) + : AdditionalData(selected) {} @@ -433,10 +412,8 @@ SparseBlockVanka::SparseBlockVanka( const SparseMatrix &M, const std::vector & selected, const unsigned int n_blocks, - const BlockingStrategy blocking_strategy, - const bool conserve_memory, - const unsigned int n_threads) - : SparseVanka(M, selected, conserve_memory, n_threads) + const BlockingStrategy blocking_strategy) + : SparseVanka(M, selected) , n_blocks(n_blocks) , dof_masks(n_blocks, std::vector(M.m(), false)) { @@ -444,6 +421,18 @@ SparseBlockVanka::SparseBlockVanka( } +template +SparseBlockVanka::SparseBlockVanka( + const SparseMatrix &M, + const std::vector & selected, + const unsigned int n_blocks, + const BlockingStrategy blocking_strategy, + const bool /*conserve_memory*/, + const unsigned int /*n_threads*/) + : SparseBlockVanka(M, selected, n_blocks, blocking_strategy) +{} + + template void SparseBlockVanka::compute_dof_masks( -- 2.39.5