Loading [MathJax]/extensions/TeX/newcommand.js
 deal.II version GIT relicensing-3067-g24f3489fdc 2025-04-14 21:10:00+00:00
\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\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
mg_coarse.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2002 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_mg_coarse_h
16#define dealii_mg_coarse_h
17
18
19#include <deal.II/base/config.h>
20
22
26
28
30
40template <typename VectorType = Vector<double>>
42{
43public:
48
53
57 void
59
63 void
65
69 void
70 operator()(const unsigned int level,
71 VectorType &dst,
72 const VectorType &src) const override;
73
74private:
81};
82
83
84
89template <class VectorType, class MatrixType>
91{
92public:
97
102
106 void
108
112 void
113 initialize(const MatrixType &matrix);
114
118 virtual void
119 operator()(const unsigned int level,
120 VectorType &dst,
121 const VectorType &src) const override;
122
123private:
127 ObserverPointer<const MatrixType,
130};
131
132
133
140template <typename VectorType,
141 typename SolverType,
142 typename MatrixType,
143 typename PreconditionerType>
145{
146public:
151
157 const MatrixType &matrix,
158 const PreconditionerType &precondition);
159
164 void
165 initialize(SolverType &solver,
166 const MatrixType &matrix,
167 const PreconditionerType &precondition);
168
172 void
174
179 virtual void
180 operator()(const unsigned int level,
181 VectorType &dst,
182 const VectorType &src) const override;
183
184private:
188 ObserverPointer<SolverType,
190 SolverType,
191 MatrixType,
192 PreconditionerType>>
194
198 ObserverPointer<const MatrixType,
200 SolverType,
201 MatrixType,
202 PreconditionerType>>
204
208 ObserverPointer<const PreconditionerType,
210 SolverType,
211 MatrixType,
212 PreconditionerType>>
214};
215
216
217
226template <typename number = double, typename VectorType = Vector<number>>
228{
229public:
234
238 void
240
241 void
242 operator()(const unsigned int level,
243 VectorType &dst,
244 const VectorType &src) const override;
245
246private:
251};
252
259template <typename number = double, typename VectorType = Vector<number>>
260class MGCoarseGridSVD : public MGCoarseGridBase<VectorType>
261{
262public:
266 MGCoarseGridSVD() = default;
267
271 void
272 initialize(const FullMatrix<number> &A, const double threshold = 0);
273
274 void
275 operator()(const unsigned int level,
276 VectorType &dst,
277 const VectorType &src) const;
278
282 void
283 log() const;
284
285private:
290};
291
294#ifndef DOXYGEN
295/* ------------------ Functions for MGCoarseGridApplySmoother -----------*/
296template <typename VectorType>
298 : coarse_smooth(nullptr)
299{}
300
301template <typename VectorType>
303 const MGSmootherBase<VectorType> &coarse_smooth)
304 : coarse_smooth(nullptr)
305{
306 initialize(coarse_smooth);
307}
308
309
310template <typename VectorType>
311void
313 const MGSmootherBase<VectorType> &coarse_smooth_)
314{
317 &coarse_smooth_, typeid(*this).name());
318}
319
320
321template <typename VectorType>
322void
324{
325 coarse_smooth = nullptr;
326}
327
328
329template <typename VectorType>
330void
332 VectorType &dst,
333 const VectorType &src) const
334{
335 coarse_smooth->apply(level, dst, src);
336}
337
338/* ------------------ Functions for MGCoarseGridApplyOperator ----------*/
339
340template <class VectorType, class PreconditionerType>
342 PreconditionerType>::MGCoarseGridApplyOperator()
343 : matrix(0, typeid(*this).name())
344{}
345
346
347
348template <class VectorType, class PreconditionerType>
350 MGCoarseGridApplyOperator(const PreconditionerType &matrix)
351 : matrix(&matrix, typeid(*this).name())
352{}
353
354
355
356template <class VectorType, class PreconditionerType>
357void
359 const PreconditionerType &matrix_)
360{
361 matrix = &matrix_;
362}
363
364
365
366template <class VectorType, class PreconditionerType>
367void
369{
370 matrix = 0;
371}
372
373
374template <class VectorType, class PreconditionerType>
375void
377 const unsigned int /*level*/,
378 VectorType &dst,
379 const VectorType &src) const
380{
381 Assert(matrix != nullptr, ExcNotInitialized());
382
383 matrix->vmult(dst, src);
384}
385
386/* ------------------ Functions for MGCoarseGridIterativeSolver ------------ */
387
388template <typename VectorType,
389 typename SolverType,
390 typename MatrixType,
391 typename PreconditionerType>
393 SolverType,
395 PreconditionerType>::MGCoarseGridIterativeSolver()
396 : solver(0, typeid(*this).name())
397 , matrix(0, typeid(*this).name())
398 , preconditioner(0, typeid(*this).name())
399{}
400
401
402
403template <typename VectorType,
404 typename SolverType,
405 typename MatrixType,
406 typename PreconditionerType>
408 SolverType,
410 PreconditionerType>::
411 MGCoarseGridIterativeSolver(SolverType &solver,
412 const MatrixType &matrix,
413 const PreconditionerType &preconditioner)
414 : solver(&solver, typeid(*this).name())
415 , matrix(&matrix, typeid(*this).name())
416 , preconditioner(&preconditioner, typeid(*this).name())
417{}
418
419
420
421template <typename VectorType,
422 typename SolverType,
423 typename MatrixType,
424 typename PreconditionerType>
425void
428 SolverType,
430 PreconditionerType>::initialize(SolverType &solver_,
431 const MatrixType &matrix_,
432 const PreconditionerType &preconditioner_)
433{
434 solver = &solver_;
435 matrix = &matrix_;
436 preconditioner = &preconditioner_;
437}
438
439
440
441template <typename VectorType,
442 typename SolverType,
443 typename MatrixType,
444 typename PreconditionerType>
445void
447 SolverType,
449 PreconditionerType>::clear()
450{
451 solver = 0;
452 matrix = 0;
453 preconditioner = 0;
454}
455
456
457
458template <typename VectorType,
459 typename SolverType,
460 typename MatrixType,
461 typename PreconditionerType>
462void
465 SolverType,
467 PreconditionerType>::operator()(const unsigned int /*level*/,
468 VectorType &dst,
469 const VectorType &src) const
470{
471 Assert(solver != nullptr, ExcNotInitialized());
472 Assert(matrix != nullptr, ExcNotInitialized());
473 Assert(preconditioner != nullptr, ExcNotInitialized());
474
475 dst = 0;
476
477 if constexpr (std::is_same_v<VectorType, typename SolverType::vector_type>)
478 {
479 solver->solve(*matrix, dst, src, *preconditioner);
480 }
481 else
482 {
483 typename SolverType::vector_type src_;
484 typename SolverType::vector_type dst_;
485
486 src_ = src;
487 dst_ = dst;
488
489 solver->solve(*matrix, dst_, src_, *preconditioner);
490
491 dst = dst_;
492 }
493}
494
495
496
497/* ------------------ Functions for MGCoarseGridHouseholder ------------ */
498
499template <typename number, typename VectorType>
501 const FullMatrix<number> *A)
502{
503 if (A != nullptr)
504 householder.initialize(*A);
505}
506
507
508
509template <typename number, typename VectorType>
510void
512 const FullMatrix<number> &A)
513{
514 householder.initialize(A);
515}
516
517
518
519template <typename number, typename VectorType>
520void
522 const unsigned int /*level*/,
523 VectorType &dst,
524 const VectorType &src) const
525{
526 householder.least_squares(dst, src);
527}
528
529//---------------------------------------------------------------------------
530
531
532
533template <typename number, typename VectorType>
534void
536 double threshold)
537{
538 matrix.reinit(A.n_rows(), A.n_cols());
539 matrix = A;
540 matrix.compute_inverse_svd(threshold);
541}
542
543
544template <typename number, typename VectorType>
545void
546MGCoarseGridSVD<number, VectorType>::operator()(const unsigned int /*level*/,
547 VectorType &dst,
548 const VectorType &src) const
549{
550 matrix.vmult(dst, src);
551}
552
553
554template <typename number, typename VectorType>
555void
557{
558 const unsigned int n = std::min(matrix.n_rows(), matrix.n_cols());
559
560 for (unsigned int i = 0; i < n; ++i)
561 deallog << ' ' << matrix.singular_value(i);
562 deallog << std::endl;
563}
564
565
566#endif // DOXYGEN
567
569
570#endif
void initialize(const MatrixType &matrix)
virtual void operator()(const unsigned int level, VectorType &dst, const VectorType &src) const override
ObserverPointer< const MatrixType, MGCoarseGridApplyOperator< VectorType, MatrixType > > matrix
Definition mg_coarse.h:129
MGCoarseGridApplyOperator(const MatrixType &matrix)
void initialize(const MGSmootherBase< VectorType > &coarse_smooth)
void operator()(const unsigned int level, VectorType &dst, const VectorType &src) const override
MGCoarseGridApplySmoother(const MGSmootherBase< VectorType > &coarse_smooth)
ObserverPointer< const MGSmootherBase< VectorType >, MGCoarseGridApplySmoother< VectorType > > coarse_smooth
Definition mg_coarse.h:80
void initialize(const FullMatrix< number > &A)
Householder< number > householder
Definition mg_coarse.h:250
MGCoarseGridHouseholder(const FullMatrix< number > *A=nullptr)
void operator()(const unsigned int level, VectorType &dst, const VectorType &src) const override
ObserverPointer< const MatrixType, MGCoarseGridIterativeSolver< VectorType, SolverType, MatrixType, PreconditionerType > > matrix
Definition mg_coarse.h:203
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)
ObserverPointer< SolverType, MGCoarseGridIterativeSolver< VectorType, SolverType, MatrixType, PreconditionerType > > solver
Definition mg_coarse.h:193
ObserverPointer< const PreconditionerType, MGCoarseGridIterativeSolver< VectorType, SolverType, MatrixType, PreconditionerType > > preconditioner
Definition mg_coarse.h:213
void log() const
MGCoarseGridSVD()=default
LAPACKFullMatrix< number > matrix
Definition mg_coarse.h:289
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:35
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
unsigned int level
Definition grid_out.cc:4635
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotInitialized()
LogStream deallog
Definition logstream.cc:38
@ matrix
Contents is actually a matrix.
constexpr char A
Tpetra::Vector< Number, LO, GO, NodeType< MemorySpace > > VectorType
Tpetra::CrsMatrix< Number, LO, GO, NodeType< MemorySpace > > MatrixType
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)