Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
sparse_vanka.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1999 - 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_sparse_vanka_h
17 #define dealii_sparse_vanka_h
18 
19 
20 
21 #include <deal.II/base/config.h>
22 
23 #include <deal.II/base/multithread_info.h>
24 #include <deal.II/base/smartpointer.h>
25 
26 #include <map>
27 #include <vector>
28 
29 DEAL_II_NAMESPACE_OPEN
30 
31 template <typename number>
32 class FullMatrix;
33 template <typename number>
34 class SparseMatrix;
35 template <typename number>
36 class Vector;
37 
38 template <typename number>
40 template <typename number>
42 
136 template <typename number>
137 class SparseVanka
138 {
139 public:
144 
151  SparseVanka();
152 
179  const std::vector<bool> & selected,
180  const bool conserve_memory = false,
181  const unsigned int n_threads = MultithreadInfo::n_threads());
182 
186  ~SparseVanka();
187 
192  {
193  public:
197  AdditionalData(const std::vector<bool> &selected,
198  const bool conserve_memory = false,
199  const unsigned int n_threads = MultithreadInfo::n_threads());
200 
204  const std::vector<bool> &selected;
205 
209  const bool conserve_mem;
210 
215  const unsigned int n_threads;
216  };
217 
218 
229  void
231  const AdditionalData & additional_data);
232 
237  template <typename number2>
238  void
239  vmult(Vector<number2> &dst, const Vector<number2> &src) const;
240 
245  template <typename number2>
246  void
247  Tvmult(Vector<number2> &dst, const Vector<number2> &src) const;
248 
256  size_type
257  m() const;
258 
266  size_type
267  n() const;
268 
269 protected:
291  template <typename number2>
292  void
293  apply_preconditioner(Vector<number2> & dst,
294  const Vector<number2> & src,
295  const std::vector<bool> *const dof_mask = nullptr) const;
296 
301  std::size_t
302  memory_consumption() const;
303 
304 private:
309 
314 
318  const std::vector<bool> *selected;
319 
324  unsigned int n_threads;
325 
330  mutable std::vector<SmartPointer<FullMatrix<float>, SparseVanka<number>>>
332 
337 
342 
346  void
348 
355  void
356  compute_inverses(const size_type begin, const size_type end);
357 
365  void
366  compute_inverse(const size_type row, std::vector<size_type> &local_indices);
367 
380  template <typename T>
381  friend class SparseBlockVanka;
382 };
383 
384 
385 
519 template <typename number>
520 class SparseBlockVanka : public SparseVanka<number>
521 {
522 public:
527 
533  {
542  };
543 
548  const std::vector<bool> & selected,
549  const unsigned int n_blocks,
550  const BlockingStrategy blocking_strategy,
551  const bool conserve_memory = false,
552  const unsigned int n_threads = MultithreadInfo::n_threads());
553 
557  template <typename number2>
558  void
559  vmult(Vector<number2> &dst, const Vector<number2> &src) const;
560 
565  std::size_t
566  memory_consumption() const;
567 
568 private:
572  const unsigned int n_blocks;
573 
581  std::vector<std::vector<bool>> dof_masks;
582 
587  void
589  const std::vector<bool> & selected,
590  const BlockingStrategy blocking_strategy);
591 };
592 
594 /* ---------------------------------- Inline functions ------------------- */
595 
596 #ifndef DOXYGEN
597 
598 template <typename number>
599 inline typename SparseVanka<number>::size_type
601 {
602  Assert(_m != 0, ExcNotInitialized());
603  return _m;
604 }
605 
606 template <typename number>
607 inline typename SparseVanka<number>::size_type
609 {
610  Assert(_n != 0, ExcNotInitialized());
611  return _n;
612 }
613 
614 template <typename number>
615 template <typename number2>
616 inline void
617 SparseVanka<number>::Tvmult(Vector<number2> & /*dst*/,
618  const Vector<number2> & /*src*/) const
619 {
620  AssertThrow(false, ExcNotImplemented());
621 }
622 
623 #endif // DOXYGEN
624 
625 DEAL_II_NAMESPACE_CLOSE
626 
627 #endif
void vmult(Vector< number2 > &dst, const Vector< number2 > &src) const
bool conserve_mem
Definition: sparse_vanka.h:313
std::size_t memory_consumption() const
void vmult(Vector< number2 > &dst, const Vector< number2 > &src) const
const unsigned int n_blocks
Definition: sparse_vanka.h:572
SparseBlockVanka(const SparseMatrix< number > &M, const std::vector< bool > &selected, const unsigned int n_blocks, const BlockingStrategy blocking_strategy, const bool conserve_memory=false, const unsigned int n_threads=MultithreadInfo::n_threads())
std::size_t memory_consumption() const
static ::ExceptionBase & ExcNotInitialized()
size_type n() const
#define AssertThrow(cond, exc)
Definition: exceptions.h:1519
void initialize(const SparseMatrix< number > &M, const AdditionalData &additional_data)
LinearAlgebra::distributed::Vector< Number > Vector
std::vector< std::vector< bool > > dof_masks
Definition: sparse_vanka.h:581
void compute_dof_masks(const SparseMatrix< number > &M, const std::vector< bool > &selected, const BlockingStrategy blocking_strategy)
types::global_dof_index size_type
Definition: sparse_vanka.h:143
#define Assert(cond, exc)
Definition: exceptions.h:1407
const std::vector< bool > * selected
Definition: sparse_vanka.h:318
unsigned int n_threads
Definition: sparse_vanka.h:324
void compute_inverses()
AdditionalData(const std::vector< bool > &selected, const bool conserve_memory=false, const unsigned int n_threads=MultithreadInfo::n_threads())
size_type m() const
const unsigned int n_threads
Definition: sparse_vanka.h:215
unsigned int global_dof_index
Definition: types.h:89
std::vector< SmartPointer< FullMatrix< float >, SparseVanka< number > > > inverses
Definition: sparse_vanka.h:331
size_type _n
Definition: sparse_vanka.h:341
void compute_inverse(const size_type row, std::vector< size_type > &local_indices)
static ::ExceptionBase & ExcNotImplemented()
static unsigned int n_threads()
SmartPointer< const SparseMatrix< number >, SparseVanka< number > > matrix
Definition: sparse_vanka.h:308
const std::vector< bool > & selected
Definition: sparse_vanka.h:204
size_type _m
Definition: sparse_vanka.h:336
void apply_preconditioner(Vector< number2 > &dst, const Vector< number2 > &src, const std::vector< bool > *const dof_mask=nullptr) const
void Tvmult(Vector< number2 > &dst, const Vector< number2 > &src) const