Reference documentation for deal.II version Git 7026f387cc 2019-10-15 14:19:01 -0400
\(\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 // Forward declarations
32 #ifndef DOXYGEN
33 template <typename number>
34 class FullMatrix;
35 template <typename number>
36 class SparseMatrix;
37 template <typename number>
38 class Vector;
39 
40 template <typename number>
41 class SparseVanka;
42 template <typename number>
43 class SparseBlockVanka;
44 #endif
45 
139 template <typename number>
141 {
142 public:
147 
154  SparseVanka();
155 
182  const std::vector<bool> & selected,
183  const bool conserve_memory = false,
184  const unsigned int n_threads = MultithreadInfo::n_threads());
185 
189  ~SparseVanka();
190 
195  {
196  public:
200  AdditionalData(const std::vector<bool> &selected,
201  const bool conserve_memory = false,
202  const unsigned int n_threads = MultithreadInfo::n_threads());
203 
207  const std::vector<bool> &selected;
208 
212  const bool conserve_mem;
213 
218  const unsigned int n_threads;
219  };
220 
221 
232  void
234  const AdditionalData & additional_data);
235 
240  template <typename number2>
241  void
242  vmult(Vector<number2> &dst, const Vector<number2> &src) const;
243 
248  template <typename number2>
249  void
250  Tvmult(Vector<number2> &dst, const Vector<number2> &src) const;
251 
259  size_type
260  m() const;
261 
269  size_type
270  n() const;
271 
272 protected:
294  template <typename number2>
295  void
296  apply_preconditioner(Vector<number2> & dst,
297  const Vector<number2> & src,
298  const std::vector<bool> *const dof_mask = nullptr) const;
299 
304  std::size_t
305  memory_consumption() const;
306 
307 private:
312 
317 
321  const std::vector<bool> *selected;
322 
327  unsigned int n_threads;
328 
333  mutable std::vector<SmartPointer<FullMatrix<float>, SparseVanka<number>>>
335 
340 
345 
349  void
351 
358  void
359  compute_inverses(const size_type begin, const size_type end);
360 
368  void
369  compute_inverse(const size_type row, std::vector<size_type> &local_indices);
370 
371  // Make the derived class a friend. This seems silly, but is actually
372  // necessary, since derived classes can only access non-public members
373  // through their @p this pointer, but not access these members as member
374  // functions of other objects of the type of this base class (i.e. like
375  // <tt>x.f()</tt>, where @p x is an object of the base class, and @p f one
376  // of it's non-public member functions).
377  //
378  // Now, in the case of the @p SparseBlockVanka class, we would like to take
379  // the address of a function of the base class in order to call it through
380  // the multithreading framework, so the derived class has to be a friend.
381  template <typename T>
382  friend class SparseBlockVanka;
383 };
384 
385 
386 
520 template <typename number>
521 class SparseBlockVanka : public SparseVanka<number>
522 {
523 public:
528 
534  {
542  adaptive
543  };
544 
549  const std::vector<bool> & selected,
550  const unsigned int n_blocks,
551  const BlockingStrategy blocking_strategy,
552  const bool conserve_memory = false,
553  const unsigned int n_threads = MultithreadInfo::n_threads());
554 
558  template <typename number2>
559  void
560  vmult(Vector<number2> &dst, const Vector<number2> &src) const;
561 
566  std::size_t
567  memory_consumption() const;
568 
569 private:
573  const unsigned int n_blocks;
574 
582  std::vector<std::vector<bool>> dof_masks;
583 
588  void
589  compute_dof_masks(const SparseMatrix<number> &M,
590  const std::vector<bool> & selected,
591  const BlockingStrategy blocking_strategy);
592 };
593 
595 /* ---------------------------------- Inline functions ------------------- */
596 
597 #ifndef DOXYGEN
598 
599 template <typename number>
600 inline typename SparseVanka<number>::size_type
602 {
603  Assert(_m != 0, ExcNotInitialized());
604  return _m;
605 }
606 
607 template <typename number>
608 inline typename SparseVanka<number>::size_type
610 {
611  Assert(_n != 0, ExcNotInitialized());
612  return _n;
613 }
614 
615 template <typename number>
616 template <typename number2>
617 inline void
618 SparseVanka<number>::Tvmult(Vector<number2> & /*dst*/,
619  const Vector<number2> & /*src*/) const
620 {
621  AssertThrow(false, ExcNotImplemented());
622 }
623 
624 #endif // DOXYGEN
625 
626 DEAL_II_NAMESPACE_CLOSE
627 
628 #endif
bool conserve_mem
Definition: sparse_vanka.h:316
void vmult(Vector< number2 > &dst, const Vector< number2 > &src) const
const unsigned int n_blocks
Definition: sparse_vanka.h:573
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:582
types::global_dof_index size_type
Definition: sparse_vanka.h:146
#define Assert(cond, exc)
Definition: exceptions.h:1407
const std::vector< bool > * selected
Definition: sparse_vanka.h:321
unsigned int n_threads
Definition: sparse_vanka.h:327
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:218
unsigned int global_dof_index
Definition: types.h:89
std::vector< SmartPointer< FullMatrix< float >, SparseVanka< number > > > inverses
Definition: sparse_vanka.h:334
size_type _n
Definition: sparse_vanka.h:344
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:311
const std::vector< bool > & selected
Definition: sparse_vanka.h:207
size_type _m
Definition: sparse_vanka.h:339
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