Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
block_matrix.h
1 // ---------------------------------------------------------------------
2 // @f$Id: block_matrix.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 2000 - 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__block_matrix_h
18 #define __deal2__block_matrix_h
19 
20 
21 #include <deal.II/base/config.h>
23 #include <deal.II/base/smartpointer.h>
24 #include <deal.II/lac/block_vector.h>
25 
26 DEAL_II_NAMESPACE_OPEN
27 
46 template <class MATRIX>
48 {
49 public:
55  BlockDiagonalMatrix (const MATRIX &M,
56  const unsigned int n_blocks);
57 
61  template <typename number1, typename number2>
62  void vmult (BlockVector<number1> &dst,
63  const BlockVector<number2> &src) const;
64 
68  template <typename number1, typename number2>
69  void Tvmult (BlockVector<number1> &dst,
70  const BlockVector<number2> &src) const;
71 private:
75  unsigned int num_blocks;
76 
81 };
82 
84 //---------------------------------------------------------------------------
85 
86 template <class MATRIX>
88  const unsigned int num_blocks)
89  :
90  num_blocks (num_blocks),
91  matrix(&M)
92 {}
93 
94 
95 template <class MATRIX>
96 template <typename number1, typename number2>
97 void
99  const BlockVector<number2> &src) const
100 {
101  Assert (dst.n_blocks()==num_blocks,
102  ExcDimensionMismatch(dst.n_blocks(),num_blocks));
103  Assert (src.n_blocks()==num_blocks,
104  ExcDimensionMismatch(src.n_blocks(),num_blocks));
105 
106  for (unsigned int i=0; i<num_blocks; ++i)
107  matrix->vmult (dst.block(i), src.block(i));
108 }
109 
110 
111 template <class MATRIX>
112 template <typename number1, typename number2>
113 void
115  const BlockVector<number2> &src) const
116 {
117  Assert (dst.n_blocks()==num_blocks,
118  ExcDimensionMismatch(dst.n_blocks(),num_blocks));
119  Assert (src.n_blocks()==num_blocks,
120  ExcDimensionMismatch(src.n_blocks(),num_blocks));
121 
122  for (unsigned int i=0; i<num_blocks; ++i)
123  matrix->Tvmult (dst.block(i), src.block(i));
124 }
125 
126 
127 DEAL_II_NAMESPACE_CLOSE
128 
129 #endif
BlockDiagonalMatrix(const MATRIX &M, const unsigned int n_blocks)
Definition: block_matrix.h:87
void Tvmult(BlockVector< number1 > &dst, const BlockVector< number2 > &src) const
Definition: block_matrix.h:114
void vmult(BlockVector< number1 > &dst, const BlockVector< number2 > &src) const
Definition: block_matrix.h:98
#define Assert(cond, exc)
Definition: exceptions.h:299
BlockType & block(const unsigned int i)
SmartPointer< const MATRIX, BlockDiagonalMatrix< MATRIX > > matrix
Definition: block_matrix.h:80
unsigned int num_blocks
Definition: block_matrix.h:75
::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)