Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
matrix_lib.templates.h
1 // ---------------------------------------------------------------------
2 // @f$Id: matrix_lib.templates.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 2002 - 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__matrix_lib_templates_h
18 #define __deal2__matrix_lib_templates_h
19 
20 #include <deal.II/lac/matrix_lib.h>
21 #include <deal.II/lac/vector.h>
22 #include <deal.II/lac/block_vector.h>
23 
25 
26 template <typename number>
27 void
29 {
30  number mean = v.mean_value();
31 
32  for (size_type i=0; i<v.size(); ++i)
33  v(i) -= mean;
34 }
35 
36 
37 
38 template <typename number>
39 void
41  const Vector<number> &src) const
42 {
43  Assert (dst.size() == src.size(),
44  ExcDimensionMismatch(dst.size(), src.size()));
45 
46  number mean = src.mean_value();
47 
48  for (size_type i=0; i<dst.size(); ++i)
49  dst(i) = src(i) - mean;
50 }
51 
52 
53 
54 template <typename number>
55 void
57  const Vector<number> &src) const
58 {
59  Assert (dst.size() == src.size(),
60  ExcDimensionMismatch(dst.size(), src.size()));
61 
62  number mean = src.mean_value();
63 
64  for (size_type i=0; i<dst.size(); ++i)
65  dst(i) += src(i) - mean;
66 }
67 
68 
69 
70 template <typename number>
71 void
73 {
76 
77  for (unsigned int i=0; i<v.n_blocks(); ++i)
78  if (i == component)
79  vmult(v.block(i), v.block(i));
80  else
81  v.block(i) = v.block(i);
82 }
83 
84 
85 
86 template <typename number>
87 void
89  const BlockVector<number> &src) const
90 {
93 
94  Assert (dst.n_blocks() == src.n_blocks(),
95  ExcDimensionMismatch(dst.n_blocks(), src.n_blocks()));
96 
97  for (unsigned int i=0; i<dst.n_blocks(); ++i)
98  if (i == component)
99  vmult(dst.block(i), src.block(i));
100  else
101  dst.block(i) = src.block(i);
102 }
103 
104 
105 
106 template <typename number>
107 void
109  const BlockVector<number> &src) const
110 {
113 
114  Assert (dst.n_blocks() == src.n_blocks(),
115  ExcDimensionMismatch(dst.n_blocks(), src.n_blocks()));
116 
117  for (unsigned int i=0; i<dst.n_blocks(); ++i)
118  if (i == component)
119  vmult_add(dst.block(i), src.block(i));
120  else
121  dst.block(i).add(src.block(i));
122 }
123 
124 
125 //----------------------------------------------------------------------//
126 
127 
128 template <class VECTOR>
130  SolverControl &c,
132  :
133  mem(m),
134  solver(c,m),
135  matrix(0),
136  precondition(0)
137 {}
138 
139 
140 template <class VECTOR>
142 {
143  if (matrix != 0) delete matrix;
144  if (precondition != 0) delete precondition;
145 }
146 
147 
148 template <class VECTOR>
149 void
150 InverseMatrixRichardson<VECTOR>::vmult(VECTOR &dst, const VECTOR &src) const
151 {
152  Assert (matrix != 0, ExcNotInitialized());
153  Assert (precondition != 0, ExcNotInitialized());
154  dst = 0.;
155  try
156  {
157  solver.solve(*matrix, dst, src, *precondition);
158  }
159  catch (...)
160  {}
161 }
162 
163 
164 
165 template <class VECTOR>
166 void
167 InverseMatrixRichardson<VECTOR>::vmult_add(VECTOR &dst, const VECTOR &src) const
168 {
169  Assert (matrix != 0, ExcNotInitialized());
170  Assert (precondition != 0, ExcNotInitialized());
171  VECTOR *aux = mem.alloc();
172  aux->reinit(dst);
173  try
174  {
175  solver.solve(*matrix, *aux, src, *precondition);
176  }
177  catch (...)
178  {}
179  dst += *aux;
180  mem.free(aux);
181 }
182 
183 
184 
185 template <class VECTOR>
186 void
187 InverseMatrixRichardson<VECTOR>::Tvmult(VECTOR &dst, const VECTOR &src) const
188 {
189  Assert (matrix != 0, ExcNotInitialized());
190  Assert (precondition != 0, ExcNotInitialized());
191  dst = 0.;
192  try
193  {
194  solver.Tsolve(*matrix, dst, src, *precondition);
195  }
196  catch (...)
197  {}
198 }
199 
200 
201 
202 template <class VECTOR>
203 void
204 InverseMatrixRichardson<VECTOR>::Tvmult_add(VECTOR &dst, const VECTOR &src) const
205 {
206  Assert (matrix != 0, ExcNotInitialized());
207  Assert (precondition != 0, ExcNotInitialized());
208  VECTOR *aux = mem.alloc();
209  aux->reinit(dst);
210  try
211  {
212  solver.Tsolve(*matrix, *aux, src, *precondition);
213  }
214  catch (...)
215  {}
216  dst += *aux;
217  mem.free(aux);
218 }
219 
220 
221 
222 
223 DEAL_II_NAMESPACE_CLOSE
224 
225 #endif
void vmult(VECTOR &, const VECTOR &) const
static const unsigned int invalid_unsigned_int
Definition: types.h:191
void vmult(Vector< number > &dst, const Vector< number > &src) const
void filter(Vector< number > &v) const
#define Assert(cond, exc)
Definition: exceptions.h:299
void Tvmult_add(VECTOR &, const VECTOR &) const
BlockType & block(const unsigned int i)
InverseMatrixRichardson(SolverControl &control, VectorMemory< VECTOR > &mem)
types::global_dof_index size_type
Definition: matrix_lib.h:337
void vmult_add(VECTOR &, const VECTOR &) const
unsigned int n_blocks() const
::ExceptionBase & ExcNotInitialized()
::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
void vmult_add(Vector< number > &dst, const Vector< number > &src) const
size_type component
Definition: matrix_lib.h:408
Number mean_value() const
std::size_t size() const
void Tvmult(VECTOR &, const VECTOR &) const