Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
sparse_vanka.h
1 // ---------------------------------------------------------------------
2 // @f$Id: sparse_vanka.h 31076 2013-10-02 17:52:30Z heister @f$
3 //
4 // Copyright (C) 1999 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__sparse_vanka_h
18 #define __deal2__sparse_vanka_h
19 
20 
21 
22 #include <deal.II/base/config.h>
23 #include <deal.II/base/smartpointer.h>
24 #include <deal.II/base/multithread_info.h>
25 
26 #include <vector>
27 #include <map>
28 
29 DEAL_II_NAMESPACE_OPEN
30 
31 template <typename number> class FullMatrix;
32 template <typename number> class SparseMatrix;
33 template <typename number> class Vector;
34 
35 template <typename number> class SparseVanka;
36 template <typename number> class SparseBlockVanka;
37 
133 template<typename number>
134 class SparseVanka
135 {
136 public:
141 
195  const std::vector<bool> &selected,
196  const bool conserve_memory = false,
197  const unsigned int n_threads = multithread_info.n_threads());
198 
203  ~SparseVanka();
204 
211  template<typename number2>
212  void vmult (Vector<number2> &dst,
213  const Vector<number2> &src) const;
214 
215 protected:
257  template<typename number2>
259  const Vector<number2> &src,
260  const std::vector<bool> *const dof_mask = 0) const;
261 
267  std::size_t memory_consumption () const;
268 
269 private:
274 
278  const bool conserve_mem;
279 
284  const std::vector<bool> &selected;
285 
292  const unsigned int n_threads;
293 
300  mutable std::vector<SmartPointer<FullMatrix<float>,SparseVanka<number> > > inverses;
301 
306  void compute_inverses ();
307 
319  void compute_inverses (const size_type begin,
320  const size_type end);
321 
336  void compute_inverse (const size_type row,
337  std::vector<size_type> &local_indices);
338 
364  template <typename T> friend class SparseBlockVanka;
365 };
366 
367 
368 
516 template<typename number>
517 class SparseBlockVanka : public SparseVanka<number>
518 {
519 public:
524 
532  {
533  index_intervals, adaptive
534  };
535 
542  const std::vector<bool> &selected,
543  const unsigned int n_blocks,
544  const BlockingStrategy blocking_strategy,
545  const bool conserve_memory = false,
546  const unsigned int n_threads = multithread_info.n_threads());
547 
551  template<typename number2>
552  void vmult (Vector<number2> &dst,
553  const Vector<number2> &src) const;
554 
560  std::size_t memory_consumption () const;
561 
562 private:
566  const unsigned int n_blocks;
567 
582  std::vector<std::vector<bool> > dof_masks;
583 
591  const std::vector<bool> &selected,
592  const BlockingStrategy blocking_strategy);
593 };
594 
596 DEAL_II_NAMESPACE_CLOSE
597 
598 #endif
types::global_dof_index size_type
Definition: sparse_vanka.h:140
const bool conserve_mem
Definition: sparse_vanka.h:278
const std::vector< bool > & selected
Definition: sparse_vanka.h:284
const unsigned int n_blocks
Definition: sparse_vanka.h:566
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=multithread_info.n_threads())
types::global_dof_index size_type
Definition: sparse_vanka.h:523
void vmult(Vector< number2 > &dst, const Vector< number2 > &src) const
void apply_preconditioner(Vector< number2 > &dst, const Vector< number2 > &src, const std::vector< bool > *const dof_mask=0) const
std::vector< SmartPointer< FullMatrix< float >, SparseVanka< number > > > inverses
Definition: sparse_vanka.h:300
std::vector< std::vector< bool > > dof_masks
Definition: sparse_vanka.h:582
const unsigned int n_threads
Definition: sparse_vanka.h:292
void compute_dof_masks(const SparseMatrix< number > &M, const std::vector< bool > &selected, const BlockingStrategy blocking_strategy)
unsigned int global_dof_index
Definition: types.h:100
std::size_t memory_consumption() const
SparseVanka(const SparseMatrix< number > &M, const std::vector< bool > &selected, const bool conserve_memory=false, const unsigned int n_threads=multithread_info.n_threads())
void compute_inverse(const size_type row, std::vector< size_type > &local_indices)
MultithreadInfo multithread_info
void vmult(Vector< number2 > &dst, const Vector< number2 > &src) const
std::size_t memory_consumption() const
SmartPointer< const SparseMatrix< number >, SparseVanka< number > > matrix
Definition: sparse_vanka.h:273