Reference documentation for deal.II version 8.1.0
mg_coarse.h
1 // ---------------------------------------------------------------------
2 // @f$Id: mg_coarse.h 30040 2013-07-18 17:06:48Z 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__mg_coarse_h
18 #define __deal2__mg_coarse_h
19 
20 
21 #include <deal.II/lac/full_matrix.h>
22 #include <deal.II/lac/matrix_lib.h>
23 #include <deal.II/lac/householder.h>
24 #include <deal.II/multigrid/mg_base.h>
25 
27 
30 
42 template<class SOLVER, class VECTOR = Vector<double> >
44 {
45 public:
50 
57  template<class MATRIX, class PRECOND>
59  const MATRIX &,
60  const PRECOND &);
61 
66 
70  template<class MATRIX, class PRECOND>
71  void initialize (SOLVER &,
72  const MATRIX &,
73  const PRECOND &);
74 
78  void clear ();
79 
87  void operator() (const unsigned int level,
88  VECTOR &dst,
89  const VECTOR &src) const;
90 
97  template <class MATRIX>
98  void set_matrix (const MATRIX &);
99 
100 private:
105 
110 
115 };
116 
117 
118 
128 template<typename number = double, class VECTOR = Vector<number> >
130 {
131 public:
137 
141  void initialize (const FullMatrix<number> &A);
142 
143  void operator() (const unsigned int level,
144  VECTOR &dst,
145  const VECTOR &src) const;
146 
147 private:
152 };
153 
162 template<typename number = double, class VECTOR = Vector<number> >
163 class MGCoarseGridSVD : public MGCoarseGridBase<VECTOR>
164 {
165 public:
170  MGCoarseGridSVD ();
171 
177  void initialize (const FullMatrix<number> &A, const double threshold = 0);
178 
179  void operator() (const unsigned int level,
180  VECTOR &dst,
181  const VECTOR &src) const;
182 
186  void log () const;
187 
188 private:
189 
194 };
195 
198 #ifndef DOXYGEN
199 /* ------------------ Functions for MGCoarseGridLACIteration ------------ */
200 
201 
202 template<class SOLVER, class VECTOR>
205  :
206  solver(0, typeid(*this).name()),
207  matrix(0),
208  precondition(0)
209 {}
210 
211 
212 template<class SOLVER, class VECTOR>
213 template<class MATRIX, class PRECOND>
216  const MATRIX &m,
217  const PRECOND &p)
218  :
219  solver(&s, typeid(*this).name())
220 {
223 }
224 
225 
226 template<class SOLVER, class VECTOR>
229 {
230  clear();
231 }
232 
233 
234 template<class SOLVER, class VECTOR>
235 template<class MATRIX, class PRECOND>
236 void
239  const MATRIX &m,
240  const PRECOND &p)
241 {
242  solver = &s;
243  if (matrix)
244  delete matrix;
246  if (precondition)
247  delete precondition;
249 }
250 
251 
252 template<class SOLVER, class VECTOR>
253 void
255 ::clear()
256 {
257  solver = 0;
258  if (matrix)
259  delete matrix;
260  matrix = 0;
261  if (precondition)
262  delete precondition;
263  precondition = 0;
264 }
265 
266 
267 template<class SOLVER, class VECTOR>
268 void
270 ::operator() (const unsigned int /* level */,
271  VECTOR &dst,
272  const VECTOR &src) const
273 {
277  solver->solve(*matrix, dst, src, *precondition);
278 }
279 
280 
281 template<class SOLVER, class VECTOR>
282 template<class MATRIX>
283 void
285 ::set_matrix(const MATRIX &m)
286 {
287  if (matrix)
288  delete matrix;
290 }
291 
292 //---------------------------------------------------------------------------
293 
294 template<typename number, class VECTOR>
296  const FullMatrix<number> *A)
297 {
298  if (A != 0) householder.initialize(*A);
299 }
300 
301 
302 
303 template<typename number, class VECTOR>
304 void
306  const FullMatrix<number> &A)
307 {
308  householder.initialize(A);
309 }
310 
311 
312 
313 template<typename number, class VECTOR>
314 void
316  const unsigned int /*level*/,
317  VECTOR &dst,
318  const VECTOR &src) const
319 {
320  householder.least_squares(dst, src);
321 }
322 
323 //---------------------------------------------------------------------------
324 
325 template<typename number, class VECTOR>
326 inline
328 {}
329 
330 
331 
332 template<typename number, class VECTOR>
333 void
335  const FullMatrix<number> &A,
336  double threshold)
337 {
338  matrix.reinit(A.n_rows(), A.n_cols());
339  matrix = A;
340  matrix.compute_inverse_svd(threshold);
341 }
342 
343 
344 template<typename number, class VECTOR>
345 void
347  const unsigned int /*level*/,
348  VECTOR &dst,
349  const VECTOR &src) const
350 {
351  matrix.vmult(dst, src);
352 }
353 
354 
355 template<typename number, class VECTOR>
356 void
358 {
359  const unsigned int n = std::min(matrix.n_rows(), matrix.n_cols());
360 
361  for (unsigned int i=0; i<n; ++i)
362  deallog << ' ' << matrix.singular_value(i);
363  deallog << std::endl;
364 }
365 
366 
367 #endif // DOXYGEN
368 
369 DEAL_II_NAMESPACE_CLOSE
370 
371 #endif
void log() const
MGCoarseGridHouseholder(const FullMatrix< number > *A=0)
PointerMatrixBase< VECTOR > * precondition
Definition: mg_coarse.h:114
void initialize(const FullMatrix< number > &A)
void initialize(const FullMatrix< number > &A, const double threshold=0)
PointerMatrixBase< VECTOR > * matrix
Definition: mg_coarse.h:109
void initialize(SOLVER &, const MATRIX &, const PRECOND &)
SmartPointer< SOLVER, MGCoarseGridLACIteration< SOLVER, VECTOR > > solver
Definition: mg_coarse.h:104
#define Assert(cond, exc)
Definition: exceptions.h:299
LAPACKFullMatrix< number > matrix
Definition: mg_coarse.h:193
void operator()(const unsigned int level, VECTOR &dst, const VECTOR &src) const
Householder< number > householder
Definition: mg_coarse.h:151
void operator()(const unsigned int level, VECTOR &dst, const VECTOR &src) const
void set_matrix(const MATRIX &)
::ExceptionBase & ExcNotInitialized()
void operator()(const unsigned int level, VECTOR &dst, const VECTOR &src) const