Reference documentation for deal.II version GIT relicensing-426-g7976cfd195 2024-04-18 21:10:01+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\}}\)
Loading...
Searching...
No Matches
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
24
26
28
38template <typename VectorType = Vector<double>>
40{
41public:
46
51
55 void
57
61 void
63
67 void
68 operator()(const unsigned int level,
69 VectorType &dst,
70 const VectorType &src) const override;
71
72private:
79};
80
81
82
89template <typename VectorType,
90 typename SolverType,
91 typename MatrixType,
92 typename PreconditionerType>
94{
95public:
100
106 const MatrixType &matrix,
107 const PreconditionerType &precondition);
108
113 void
114 initialize(SolverType &solver,
115 const MatrixType &matrix,
116 const PreconditionerType &precondition);
117
121 void
123
128 virtual void
129 operator()(const unsigned int level,
130 VectorType &dst,
131 const VectorType &src) const override;
132
133private:
137 SmartPointer<SolverType,
139 SolverType,
140 MatrixType,
141 PreconditionerType>>
143
147 SmartPointer<const MatrixType,
149 SolverType,
150 MatrixType,
151 PreconditionerType>>
153
157 SmartPointer<const PreconditionerType,
159 SolverType,
160 MatrixType,
161 PreconditionerType>>
163};
164
165
166
175template <typename number = double, typename VectorType = Vector<number>>
177{
178public:
183
187 void
189
190 void
191 operator()(const unsigned int level,
192 VectorType &dst,
193 const VectorType &src) const override;
194
195private:
200};
201
208template <typename number = double, typename VectorType = Vector<number>>
209class MGCoarseGridSVD : public MGCoarseGridBase<VectorType>
210{
211public:
215 MGCoarseGridSVD() = default;
216
220 void
221 initialize(const FullMatrix<number> &A, const double threshold = 0);
222
223 void
224 operator()(const unsigned int level,
225 VectorType &dst,
226 const VectorType &src) const;
227
231 void
232 log() const;
233
234private:
239};
240
243#ifndef DOXYGEN
244/* ------------------ Functions for MGCoarseGridApplySmoother -----------*/
245template <typename VectorType>
247 : coarse_smooth(nullptr)
248{}
249
250template <typename VectorType>
252 const MGSmootherBase<VectorType> &coarse_smooth)
253 : coarse_smooth(nullptr)
254{
255 initialize(coarse_smooth);
256}
257
258
259template <typename VectorType>
260void
262 const MGSmootherBase<VectorType> &coarse_smooth_)
263{
264 coarse_smooth =
267 typeid(*this).name());
268}
269
270
271template <typename VectorType>
272void
274{
275 coarse_smooth = nullptr;
276}
277
278
279template <typename VectorType>
280void
282 VectorType &dst,
283 const VectorType &src) const
284{
285 coarse_smooth->apply(level, dst, src);
286}
287
288/* ------------------ Functions for MGCoarseGridIterativeSolver ------------ */
289
290template <typename VectorType,
291 typename SolverType,
292 typename MatrixType,
293 typename PreconditionerType>
295 SolverType,
296 MatrixType,
297 PreconditionerType>::MGCoarseGridIterativeSolver()
298 : solver(0, typeid(*this).name())
299 , matrix(0, typeid(*this).name())
300 , preconditioner(0, typeid(*this).name())
301{}
302
303
304
305template <typename VectorType,
306 typename SolverType,
307 typename MatrixType,
308 typename PreconditionerType>
310 SolverType,
311 MatrixType,
312 PreconditionerType>::
313 MGCoarseGridIterativeSolver(SolverType &solver,
314 const MatrixType &matrix,
315 const PreconditionerType &preconditioner)
316 : solver(&solver, typeid(*this).name())
317 , matrix(&matrix, typeid(*this).name())
318 , preconditioner(&preconditioner, typeid(*this).name())
319{}
320
321
322
323template <typename VectorType,
324 typename SolverType,
325 typename MatrixType,
326 typename PreconditionerType>
327void
329 VectorType,
330 SolverType,
331 MatrixType,
332 PreconditionerType>::initialize(SolverType &solver_,
333 const MatrixType &matrix_,
334 const PreconditionerType &preconditioner_)
335{
336 solver = &solver_;
337 matrix = &matrix_;
338 preconditioner = &preconditioner_;
339}
340
341
342
343template <typename VectorType,
344 typename SolverType,
345 typename MatrixType,
346 typename PreconditionerType>
347void
349 SolverType,
350 MatrixType,
351 PreconditionerType>::clear()
352{
353 solver = 0;
354 matrix = 0;
355 preconditioner = 0;
356}
357
358
359
360namespace internal
361{
363 {
364 template <typename VectorType,
365 typename SolverType,
366 typename MatrixType,
367 typename PreconditionerType,
368 std::enable_if_t<
369 std::is_same_v<VectorType, typename SolverType::vector_type>,
370 VectorType> * = nullptr>
371 void
372 solve(SolverType &solver,
373 const MatrixType &matrix,
374 const PreconditionerType &preconditioner,
375 VectorType &dst,
376 const VectorType &src)
377 {
378 solver.solve(matrix, dst, src, preconditioner);
379 }
380
381 template <typename VectorType,
382 typename SolverType,
383 typename MatrixType,
384 typename PreconditionerType,
385 std::enable_if_t<
386 !std::is_same_v<VectorType, typename SolverType::vector_type>,
387 VectorType> * = nullptr>
388 void
389 solve(SolverType &solver,
390 const MatrixType &matrix,
391 const PreconditionerType &preconditioner,
392 VectorType &dst,
393 const VectorType &src)
394 {
395 typename SolverType::vector_type src_;
396 typename SolverType::vector_type dst_;
397
398 src_ = src;
399 dst_ = dst;
400
401 solver.solve(matrix, dst_, src_, preconditioner);
402
403 dst = dst_;
404 }
405 } // namespace MGCoarseGridIterativeSolver
406} // namespace internal
407
408
409
410template <typename VectorType,
411 typename SolverType,
412 typename MatrixType,
413 typename PreconditionerType>
414void
416 VectorType,
417 SolverType,
418 MatrixType,
419 PreconditionerType>::operator()(const unsigned int /*level*/,
420 VectorType &dst,
421 const VectorType &src) const
422{
423 Assert(solver != nullptr, ExcNotInitialized());
424 Assert(matrix != nullptr, ExcNotInitialized());
425 Assert(preconditioner != nullptr, ExcNotInitialized());
426
427 dst = 0;
428 internal::MGCoarseGridIterativeSolver::solve(
429 *solver, *matrix, *preconditioner, dst, src);
430}
431
432
433
434/* ------------------ Functions for MGCoarseGridHouseholder ------------ */
435
436template <typename number, typename VectorType>
438 const FullMatrix<number> *A)
439{
440 if (A != nullptr)
441 householder.initialize(*A);
442}
443
444
445
446template <typename number, typename VectorType>
447void
449 const FullMatrix<number> &A)
450{
451 householder.initialize(A);
452}
453
454
455
456template <typename number, typename VectorType>
457void
459 const unsigned int /*level*/,
460 VectorType &dst,
461 const VectorType &src) const
462{
463 householder.least_squares(dst, src);
464}
465
466//---------------------------------------------------------------------------
467
468
469
470template <typename number, typename VectorType>
471void
473 double threshold)
474{
475 matrix.reinit(A.n_rows(), A.n_cols());
476 matrix = A;
477 matrix.compute_inverse_svd(threshold);
478}
479
480
481template <typename number, typename VectorType>
482void
483MGCoarseGridSVD<number, VectorType>::operator()(const unsigned int /*level*/,
484 VectorType &dst,
485 const VectorType &src) const
486{
487 matrix.vmult(dst, src);
488}
489
490
491template <typename number, typename VectorType>
492void
494{
495 const unsigned int n = std::min(matrix.n_rows(), matrix.n_cols());
496
497 for (unsigned int i = 0; i < n; ++i)
498 deallog << ' ' << matrix.singular_value(i);
499 deallog << std::endl;
500}
501
502
503#endif // DOXYGEN
504
506
507#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:78
MGCoarseGridApplySmoother(const MGSmootherBase< VectorType > &coarse_smooth)
void initialize(const FullMatrix< number > &A)
Householder< number > householder
Definition mg_coarse.h:199
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:152
SmartPointer< SolverType, MGCoarseGridIterativeSolver< VectorType, SolverType, MatrixType, PreconditionerType > > solver
Definition mg_coarse.h:142
SmartPointer< const PreconditionerType, MGCoarseGridIterativeSolver< VectorType, SolverType, MatrixType, PreconditionerType > > preconditioner
Definition mg_coarse.h:162
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:238
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:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
unsigned int level
Definition grid_out.cc:4626
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotInitialized()
LogStream deallog
Definition logstream.cc:36
@ matrix
Contents is actually a matrix.
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)