Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
mg_coarse.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2002 - 2022 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_mg_coarse_h
17#define dealii_mg_coarse_h
18
19
20#include <deal.II/base/config.h>
21
25
27
29
39template <class VectorType = Vector<double>>
41{
42public:
47
52
56 void
58
62 void
64
68 void
69 operator()(const unsigned int level,
70 VectorType & dst,
71 const VectorType & src) const override;
72
73private:
80};
81
82
83
90template <class VectorType,
91 class SolverType,
92 class MatrixType,
93 class PreconditionerType>
95{
96public:
101
107 const MatrixType & matrix,
108 const PreconditionerType &precondition);
109
114 void
115 initialize(SolverType & solver,
116 const MatrixType & matrix,
117 const PreconditionerType &precondition);
118
122 void
124
129 virtual void
130 operator()(const unsigned int level,
131 VectorType & dst,
132 const VectorType & src) const override;
133
134private:
138 SmartPointer<SolverType,
140 SolverType,
141 MatrixType,
142 PreconditionerType>>
144
148 SmartPointer<const MatrixType,
150 SolverType,
151 MatrixType,
152 PreconditionerType>>
154
158 SmartPointer<const PreconditionerType,
160 SolverType,
161 MatrixType,
162 PreconditionerType>>
164};
165
166
167
176template <typename number = double, class VectorType = Vector<number>>
178{
179public:
184
188 void
190
191 void
192 operator()(const unsigned int level,
193 VectorType & dst,
194 const VectorType & src) const override;
195
196private:
201};
202
209template <typename number = double, class VectorType = Vector<number>>
210class MGCoarseGridSVD : public MGCoarseGridBase<VectorType>
211{
212public:
216 MGCoarseGridSVD() = default;
217
221 void
222 initialize(const FullMatrix<number> &A, const double threshold = 0);
223
224 void
225 operator()(const unsigned int level,
226 VectorType & dst,
227 const VectorType & src) const;
228
232 void
233 log() const;
234
235private:
240};
241
244#ifndef DOXYGEN
245/* ------------------ Functions for MGCoarseGridApplySmoother -----------*/
246template <class VectorType>
248 : coarse_smooth(nullptr)
249{}
250
251template <class VectorType>
253 const MGSmootherBase<VectorType> &coarse_smooth)
254 : coarse_smooth(nullptr)
255{
256 initialize(coarse_smooth);
257}
258
259
260template <class VectorType>
261void
263 const MGSmootherBase<VectorType> &coarse_smooth_)
264{
265 coarse_smooth =
268 typeid(*this).name());
269}
270
271
272template <class VectorType>
273void
275{
276 coarse_smooth = nullptr;
277}
278
279
280template <class VectorType>
281void
283 VectorType & dst,
284 const VectorType & src) const
285{
286 coarse_smooth->apply(level, dst, src);
287}
288
289/* ------------------ Functions for MGCoarseGridIterativeSolver ------------ */
290
291template <class VectorType,
292 class SolverType,
293 class MatrixType,
294 class PreconditionerType>
296 SolverType,
297 MatrixType,
298 PreconditionerType>::MGCoarseGridIterativeSolver()
299 : solver(0, typeid(*this).name())
300 , matrix(0, typeid(*this).name())
301 , preconditioner(0, typeid(*this).name())
302{}
303
304
305
306template <class VectorType,
307 class SolverType,
308 class MatrixType,
309 class PreconditionerType>
311 SolverType,
312 MatrixType,
313 PreconditionerType>::
314 MGCoarseGridIterativeSolver(SolverType & solver,
315 const MatrixType & matrix,
316 const PreconditionerType &preconditioner)
317 : solver(&solver, typeid(*this).name())
318 , matrix(&matrix, typeid(*this).name())
319 , preconditioner(&preconditioner, typeid(*this).name())
320{}
321
322
323
324template <class VectorType,
325 class SolverType,
326 class MatrixType,
327 class PreconditionerType>
328void
330 VectorType,
331 SolverType,
332 MatrixType,
333 PreconditionerType>::initialize(SolverType & solver_,
334 const MatrixType & matrix_,
335 const PreconditionerType &preconditioner_)
336{
337 solver = &solver_;
338 matrix = &matrix_;
339 preconditioner = &preconditioner_;
340}
341
342
343
344template <class VectorType,
345 class SolverType,
346 class MatrixType,
347 class PreconditionerType>
348void
350 SolverType,
351 MatrixType,
352 PreconditionerType>::clear()
353{
354 solver = 0;
355 matrix = 0;
356 preconditioner = 0;
357}
358
359
360
361namespace internal
362{
364 {
365 template <
366 class VectorType,
367 class SolverType,
368 class MatrixType,
369 class PreconditionerType,
370 std::enable_if_t<
371 std::is_same<VectorType, typename SolverType::vector_type>::value,
372 VectorType> * = nullptr>
373 void
374 solve(SolverType & solver,
375 const MatrixType & matrix,
376 const PreconditionerType &preconditioner,
377 VectorType & dst,
378 const VectorType & src)
379 {
380 solver.solve(matrix, dst, src, preconditioner);
381 }
382
383 template <
384 class VectorType,
385 class SolverType,
386 class MatrixType,
387 class PreconditionerType,
388 std::enable_if_t<
389 !std::is_same<VectorType, typename SolverType::vector_type>::value,
390 VectorType> * = nullptr>
391 void
392 solve(SolverType & solver,
393 const MatrixType & matrix,
394 const PreconditionerType &preconditioner,
395 VectorType & dst,
396 const VectorType & src)
397 {
398 typename SolverType::vector_type src_;
399 typename SolverType::vector_type dst_;
400
401 src_ = src;
402 dst_ = dst;
403
404 solver.solve(matrix, dst_, src_, preconditioner);
405
406 dst = dst_;
407 }
408 } // namespace MGCoarseGridIterativeSolver
409} // namespace internal
410
411
412
413template <class VectorType,
414 class SolverType,
415 class MatrixType,
416 class PreconditionerType>
417void
419 VectorType,
420 SolverType,
421 MatrixType,
422 PreconditionerType>::operator()(const unsigned int /*level*/,
423 VectorType & dst,
424 const VectorType &src) const
425{
426 Assert(solver != nullptr, ExcNotInitialized());
427 Assert(matrix != nullptr, ExcNotInitialized());
428 Assert(preconditioner != nullptr, ExcNotInitialized());
429
430 dst = 0;
431 internal::MGCoarseGridIterativeSolver::solve(
432 *solver, *matrix, *preconditioner, dst, src);
433}
434
435
436
437/* ------------------ Functions for MGCoarseGridHouseholder ------------ */
438
439template <typename number, class VectorType>
441 const FullMatrix<number> *A)
442{
443 if (A != nullptr)
444 householder.initialize(*A);
445}
446
447
448
449template <typename number, class VectorType>
450void
452 const FullMatrix<number> &A)
453{
454 householder.initialize(A);
455}
456
457
458
459template <typename number, class VectorType>
460void
462 const unsigned int /*level*/,
463 VectorType & dst,
464 const VectorType &src) const
465{
466 householder.least_squares(dst, src);
467}
468
469//---------------------------------------------------------------------------
470
471
472
473template <typename number, class VectorType>
474void
476 double threshold)
477{
478 matrix.reinit(A.n_rows(), A.n_cols());
479 matrix = A;
480 matrix.compute_inverse_svd(threshold);
481}
482
483
484template <typename number, class VectorType>
485void
486MGCoarseGridSVD<number, VectorType>::operator()(const unsigned int /*level*/,
487 VectorType & dst,
488 const VectorType &src) const
489{
490 matrix.vmult(dst, src);
491}
492
493
494template <typename number, class VectorType>
495void
497{
498 const unsigned int n = std::min(matrix.n_rows(), matrix.n_cols());
499
500 for (unsigned int i = 0; i < n; ++i)
501 deallog << ' ' << matrix.singular_value(i);
502 deallog << std::endl;
503}
504
505
506#endif // DOXYGEN
507
509
510#endif
void initialize(const MGSmootherBase< VectorType > &coarse_smooth)
void operator()(const unsigned int level, VectorType &dst, const VectorType &src) const override
SmartPointer< const MGSmootherBase< VectorType >, MGCoarseGridApplySmoother< VectorType > > coarse_smooth
Definition mg_coarse.h:79
MGCoarseGridApplySmoother(const MGSmootherBase< VectorType > &coarse_smooth)
void initialize(const FullMatrix< number > &A)
Householder< number > householder
Definition mg_coarse.h:200
MGCoarseGridHouseholder(const FullMatrix< number > *A=nullptr)
void operator()(const unsigned int level, VectorType &dst, const VectorType &src) const override
SmartPointer< const MatrixType, MGCoarseGridIterativeSolver< VectorType, SolverType, MatrixType, PreconditionerType > > matrix
Definition mg_coarse.h:153
SmartPointer< SolverType, MGCoarseGridIterativeSolver< VectorType, SolverType, MatrixType, PreconditionerType > > solver
Definition mg_coarse.h:143
SmartPointer< const PreconditionerType, MGCoarseGridIterativeSolver< VectorType, SolverType, MatrixType, PreconditionerType > > preconditioner
Definition mg_coarse.h:163
void initialize(SolverType &solver, const MatrixType &matrix, const PreconditionerType &precondition)
virtual void operator()(const unsigned int level, VectorType &dst, const VectorType &src) const override
MGCoarseGridIterativeSolver(SolverType &solver, const MatrixType &matrix, const PreconditionerType &precondition)
void log() const
MGCoarseGridSVD()=default
LAPACKFullMatrix< number > matrix
Definition mg_coarse.h:239
void operator()(const unsigned int level, VectorType &dst, const VectorType &src) const
void initialize(const FullMatrix< number > &A, const double threshold=0)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
unsigned int level
Definition grid_out.cc:4618
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotInitialized()
LogStream deallog
Definition logstream.cc:37
@ matrix
Contents is actually a matrix.
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)