15#ifndef dealii_precondition_block_base_h
16#define dealii_precondition_block_base_h
35template <
typename number>
37template <
typename number>
56template <
typename number>
147 template <
typename number2>
156 template <
typename number2>
301template <
typename number>
305 , n_diagonal_blocks(0)
306 , var_store_diagonals(store)
307 , var_same_diagonal(false)
308 , var_inverses_ready(false)
312template <
typename number>
316 if (var_inverse_full.size() != 0)
317 var_inverse_full.erase(var_inverse_full.begin(), var_inverse_full.end());
318 if (var_inverse_householder.size() != 0)
319 var_inverse_householder.erase(var_inverse_householder.begin(),
320 var_inverse_householder.end());
321 if (var_inverse_svd.size() != 0)
322 var_inverse_svd.erase(var_inverse_svd.begin(), var_inverse_svd.end());
323 if (var_diagonal.size() != 0)
324 var_diagonal.erase(var_diagonal.begin(), var_diagonal.end());
325 var_same_diagonal =
false;
326 var_inverses_ready =
false;
327 n_diagonal_blocks = 0;
330template <
typename number>
338 var_same_diagonal = compress;
339 var_inverses_ready =
false;
340 n_diagonal_blocks = n;
347 var_inverse_full.resize(1);
348 var_inverse_full[0].reinit(b, b);
351 var_inverse_householder.resize(1);
354 var_inverse_svd.resize(1);
355 var_inverse_svd[0].reinit(b, b);
361 if (store_diagonals())
363 var_diagonal.resize(1);
364 var_diagonal[0].reinit(b, b);
377 if (store_diagonals())
380 var_diagonal.swap(tmp);
388 var_inverse_full.swap(tmp);
392 var_inverse_householder.resize(n);
396 std::vector<LAPACKFullMatrix<number>> tmp(
398 var_inverse_svd.swap(tmp);
408template <
typename number>
412 return n_diagonal_blocks;
415template <
typename number>
416template <
typename number2>
422 const size_type ii = same_diagonal() ? 0U : i;
428 var_inverse_full[ii].vmult(dst, src);
432 var_inverse_householder[ii].vmult(dst, src);
436 var_inverse_svd[ii].vmult(dst, src);
444template <
typename number>
445template <
typename number2>
451 const size_type ii = same_diagonal() ? 0U : i;
457 var_inverse_full[ii].Tvmult(dst, src);
461 var_inverse_householder[ii].Tvmult(dst, src);
465 var_inverse_svd[ii].Tvmult(dst, src);
473template <
typename number>
478 return var_inverse_full[0];
481 return var_inverse_full[i];
485template <
typename number>
490 return var_inverse_householder[0];
493 return var_inverse_householder[i];
497template <
typename number>
502 return var_inverse_svd[0];
505 return var_inverse_svd[i];
509template <
typename number>
513 Assert(store_diagonals(), ExcDiagonalsNotStored());
516 return var_diagonal[0];
519 return var_diagonal[i];
523template <
typename number>
527 Assert(var_inverse_full.size() != 0, ExcInverseNotAvailable());
530 return var_inverse_full[0];
533 return var_inverse_full[i];
537template <
typename number>
541 Assert(var_inverse_householder.size() != 0, ExcInverseNotAvailable());
544 return var_inverse_householder[0];
547 return var_inverse_householder[i];
551template <
typename number>
555 Assert(var_inverse_svd.size() != 0, ExcInverseNotAvailable());
558 return var_inverse_svd[0];
561 return var_inverse_svd[i];
565template <
typename number>
569 Assert(store_diagonals(), ExcDiagonalsNotStored());
572 return var_diagonal[0];
575 return var_diagonal[i];
579template <
typename number>
583 return var_same_diagonal;
587template <
typename number>
591 return var_store_diagonals;
595template <
typename number>
599 var_inverses_ready = x;
603template <
typename number>
607 return var_inverses_ready;
611template <
typename number>
615 deallog <<
"PreconditionBlockBase: " << size() <<
" blocks; ";
617 if (inversion == svd)
619 unsigned int kermin = 100000000, kermax = 0;
620 double sigmin = 1.e300, sigmax = -1.e300;
621 double kappamin = 1.e300, kappamax = -1.e300;
627 while (k <= matrix.n_cols() &&
628 matrix.singular_value(matrix.n_cols() - k) == 0)
630 const double s0 = matrix.singular_value(0);
631 const double sm = matrix.singular_value(matrix.n_cols() - k);
632 const double co = sm / s0;
647 deallog <<
"dim ker [" << kermin <<
':' << kermax <<
"] sigma [" << sigmin
648 <<
':' << sigmax <<
"] kappa [" << kappamin <<
':' << kappamax
651 else if (inversion == householder)
654 else if (inversion == gauss_jordan)
664template <
typename number>
668 std::size_t mem =
sizeof(*this);
669 for (
size_type i = 0; i < var_inverse_full.size(); ++i)
671 for (
size_type i = 0; i < var_diagonal.size(); ++i)
bool same_diagonal() const
std::vector< FullMatrix< number > > var_diagonal
unsigned int n_diagonal_blocks
std::size_t memory_consumption() const
std::vector< LAPACKFullMatrix< number > > var_inverse_svd
void inverses_computed(bool are_they)
void reinit(unsigned int nblocks, size_type blocksize, bool compress, Inversion method=gauss_jordan)
void inverse_vmult(size_type i, Vector< number2 > &dst, const Vector< number2 > &src) const
void log_statistics() const
PreconditionBlockBase(bool store_diagonals=false, Inversion method=gauss_jordan)
Householder< number > & inverse_householder(size_type i)
std::vector< FullMatrix< number > > var_inverse_full
bool inverses_ready() const
unsigned int size() const
const FullMatrix< number > & diagonal(size_type i) const
std::vector< Householder< number > > var_inverse_householder
FullMatrix< number > & inverse(size_type i)
FullMatrix< number > & diagonal(size_type i)
bool store_diagonals() const
const Householder< number > & inverse_householder(size_type i) const
const FullMatrix< number > & inverse(size_type i) const
void inverse_Tvmult(size_type i, Vector< number2 > &dst, const Vector< number2 > &src) const
LAPACKFullMatrix< number > & inverse_svd(size_type i)
const LAPACKFullMatrix< number > & inverse_svd(size_type i) const
~PreconditionBlockBase()=default
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define DeclException0(Exception0)
static ::ExceptionBase & ExcDiagonalsNotStored()
#define Assert(cond, exc)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInverseNotAvailable()
#define DEAL_II_NOT_IMPLEMENTED()
std::enable_if_t< std::is_fundamental_v< T >, std::size_t > memory_consumption(const T &t)
unsigned int global_dof_index