16 #ifndef dealii_precondition_h 17 #define dealii_precondition_h 40 template <
typename number>
42 template <
typename number>
48 template <
typename,
typename>
108 template <
typename MatrixType>
110 initialize(
const MatrixType &
matrix,
116 template <
class VectorType>
124 template <
class VectorType>
131 template <
class VectorType>
139 template <
class VectorType>
238 template <
typename MatrixType>
245 template <
class VectorType>
253 template <
class VectorType>
259 template <
class VectorType>
267 template <
class VectorType>
357 template <
typename MatrixType = SparseMatrix<
double>,
358 class VectorType = Vector<
double>>
401 template <
typename MatrixType = SparseMatrix<
double>>
433 initialize(
const MatrixType &
A,
496 template <
typename MatrixType = SparseMatrix<
double>>
509 template <
class VectorType>
517 template <
class VectorType>
524 template <
class VectorType>
531 template <
class VectorType>
582 template <
typename MatrixType = SparseMatrix<
double>>
595 template <
class VectorType>
602 template <
class VectorType>
609 template <
class VectorType>
616 template <
class VectorType>
649 template <
typename MatrixType = SparseMatrix<
double>>
676 initialize(
const MatrixType &
A,
683 template <
class VectorType>
691 template <
class VectorType>
699 template <
class VectorType>
706 template <
class VectorType>
748 template <
typename MatrixType = SparseMatrix<
double>>
774 const std::vector<size_type> &permutation,
775 const std::vector<size_type> &inverse_permutation,
806 initialize(
const MatrixType &
A,
807 const std::vector<size_type> &permutation,
808 const std::vector<size_type> &inverse_permutation,
824 initialize(
const MatrixType &A,
const AdditionalData &additional_data);
829 template <
class VectorType>
836 template <
class VectorType>
988 template <
typename MatrixType = SparseMatrix<
double>,
989 typename VectorType = Vector<
double>,
990 typename PreconditionerType = DiagonalMatrix<VectorType>>
1009 const double smoothing_range = 0.,
1010 const unsigned int eig_cg_n_iterations = 8,
1011 const double eig_cg_residual = 1
e-2,
1012 const double max_eigenvalue = 1);
1097 initialize(
const MatrixType &
matrix,
1174 , max_eigenvalue_estimate{std::numeric_limits<double>::lowest()}
1193 estimate_eigenvalues(
const VectorType &src)
const;
1261 template <
typename MatrixType>
1266 n_rows = matrix.m();
1267 n_columns = matrix.n();
1271 template <
class VectorType>
1280 template <
class VectorType>
1287 template <
class VectorType>
1296 template <
class VectorType>
1320 const double relaxation)
1321 : relaxation(relaxation)
1331 relaxation = add_data.relaxation;
1345 template <
typename MatrixType>
1348 const MatrixType & matrix,
1352 n_rows = matrix.m();
1353 n_columns = matrix.n();
1358 template <
class VectorType>
1363 std::is_same<size_type, typename VectorType::size_type>::value,
1364 "PreconditionRichardson and VectorType must have the same size_type.");
1366 dst.equ(relaxation, src);
1371 template <
class VectorType>
1376 std::is_same<size_type, typename VectorType::size_type>::value,
1377 "PreconditionRichardson and VectorType must have the same size_type.");
1379 dst.equ(relaxation, src);
1382 template <
class VectorType>
1387 std::is_same<size_type, typename VectorType::size_type>::value,
1388 "PreconditionRichardson and VectorType must have the same size_type.");
1390 dst.add(relaxation, src);
1395 template <
class VectorType>
1400 std::is_same<size_type, typename VectorType::size_type>::value,
1401 "PreconditionRichardson and VectorType must have the same size_type.");
1403 dst.add(relaxation, src);
1422 template <
typename MatrixType>
1428 relaxation = parameters.relaxation;
1432 template <
typename MatrixType>
1439 template <
typename MatrixType>
1447 template <
typename MatrixType>
1457 template <
typename MatrixType>
1458 template <
class VectorType>
1466 "PreconditionJacobi and VectorType must have the same size_type.");
1469 this->
A->precondition_Jacobi(dst, src, this->relaxation);
1474 template <
typename MatrixType>
1475 template <
class VectorType>
1483 "PreconditionJacobi and VectorType must have the same size_type.");
1486 this->
A->precondition_Jacobi(dst, src, this->relaxation);
1491 template <
typename MatrixType>
1492 template <
class VectorType>
1500 "PreconditionJacobi and VectorType must have the same size_type.");
1503 this->
A->Jacobi_step(dst, src, this->relaxation);
1508 template <
typename MatrixType>
1509 template <
class VectorType>
1517 "PreconditionJacobi and VectorType must have the same size_type.");
1526 template <
typename MatrixType>
1527 template <
class VectorType>
1533 "PreconditionSOR and VectorType must have the same size_type.");
1536 this->
A->precondition_SOR(dst, src, this->relaxation);
1541 template <
typename MatrixType>
1542 template <
class VectorType>
1549 "PreconditionSOR and VectorType must have the same size_type.");
1552 this->
A->precondition_TSOR(dst, src, this->relaxation);
1557 template <
typename MatrixType>
1558 template <
class VectorType>
1564 "PreconditionSOR and VectorType must have the same size_type.");
1567 this->
A->SOR_step(dst, src, this->relaxation);
1572 template <
typename MatrixType>
1573 template <
class VectorType>
1579 "PreconditionSOR and VectorType must have the same size_type.");
1582 this->
A->TSOR_step(dst, src, this->relaxation);
1589 template <
typename MatrixType>
1592 const MatrixType & rA,
1593 const typename BaseClass::AdditionalData ¶meters)
1607 pos_right_of_diagonal.resize(n, static_cast<std::size_t>(-1));
1615 it = mat->
begin(row) + 1;
1616 for (; it < mat->
end(row); ++it)
1617 if (it->column() > row)
1619 pos_right_of_diagonal[row] = it - mat->
begin();
1625 template <
typename MatrixType>
1626 template <
class VectorType>
1634 "PreconditionSSOR and VectorType must have the same size_type.");
1637 this->
A->precondition_SSOR(dst, src, this->relaxation, pos_right_of_diagonal);
1642 template <
typename MatrixType>
1643 template <
class VectorType>
1651 "PreconditionSSOR and VectorType must have the same size_type.");
1654 this->
A->precondition_SSOR(dst, src, this->relaxation, pos_right_of_diagonal);
1659 template <
typename MatrixType>
1660 template <
class VectorType>
1667 "PreconditionSSOR and VectorType must have the same size_type.");
1670 this->
A->SSOR_step(dst, src, this->relaxation);
1675 template <
typename MatrixType>
1676 template <
class VectorType>
1684 "PreconditionSSOR and VectorType must have the same size_type.");
1693 template <
typename MatrixType>
1696 const MatrixType & rA,
1697 const std::vector<size_type> & p,
1698 const std::vector<size_type> & ip,
1702 inverse_permutation = &ip;
1707 template <
typename MatrixType>
1713 additional_data.permutation,
1714 additional_data.inverse_permutation,
1715 additional_data.parameters);
1719 template <
typename MatrixType>
1720 template <
typename VectorType>
1728 "PreconditionPSOR and VectorType must have the same size_type.");
1732 this->A->PSOR(dst, *permutation, *inverse_permutation, this->relaxation);
1737 template <
typename MatrixType>
1738 template <
class VectorType>
1746 "PreconditionPSOR and VectorType must have the same size_type.");
1750 this->A->TPSOR(dst, *permutation, *inverse_permutation, this->relaxation);
1753 template <
typename MatrixType>
1755 const std::vector<size_type> &permutation,
1756 const std::vector<size_type> &inverse_permutation,
1758 : permutation(permutation)
1759 , inverse_permutation(inverse_permutation)
1760 , parameters(parameters)
1767 template <
typename MatrixType,
class VectorType>
1769 const MatrixType & M,
1770 const function_ptr method)
1772 , precondition(method)
1777 template <
typename MatrixType,
class VectorType>
1783 (matrix.*precondition)(dst, src);
1788 template <
typename MatrixType>
1790 const double relaxation)
1791 : relaxation(relaxation)
1800 namespace PreconditionChebyshevImplementation
1808 template <
typename VectorType,
typename PreconditionerType>
1811 const PreconditionerType &preconditioner,
1812 const unsigned int iteration_index,
1813 const double factor1,
1814 const double factor2,
1820 if (iteration_index == 0)
1822 solution.equ(factor2, rhs);
1823 preconditioner.vmult(solution_old, solution);
1825 else if (iteration_index == 1)
1828 temp_vector1.sadd(-1.0, 1.0, rhs);
1829 preconditioner.vmult(solution_old, temp_vector1);
1832 solution_old.sadd(factor2, 1 + factor1, solution);
1837 temp_vector1.sadd(-1.0, 1.0, rhs);
1838 preconditioner.vmult(temp_vector2, temp_vector1);
1841 solution_old.sadd(-factor1, factor2, temp_vector2);
1842 solution_old.add(1 + factor1, solution);
1845 solution.swap(solution_old);
1851 template <
typename Number>
1852 struct VectorUpdater
1854 VectorUpdater(
const Number * rhs,
1855 const Number * matrix_diagonal_inverse,
1856 const unsigned int iteration_index,
1857 const Number factor1,
1858 const Number factor2,
1859 Number * solution_old,
1860 Number * tmp_vector,
1863 , matrix_diagonal_inverse(matrix_diagonal_inverse)
1864 , iteration_index(iteration_index)
1868 , tmp_vector(tmp_vector)
1869 , solution(solution)
1873 apply_to_subrange(
const std::size_t
begin,
const std::size_t
end)
const 1879 const Number factor1 = this->factor1;
1880 const Number factor1_plus_1 = 1. + this->factor1;
1881 const Number factor2 = this->factor2;
1882 if (iteration_index == 0)
1885 for (std::size_t i = begin; i <
end; ++i)
1886 solution[i] = factor2 * matrix_diagonal_inverse[i] * rhs[i];
1888 else if (iteration_index == 1)
1892 for (std::size_t i = begin; i <
end; ++i)
1896 factor1_plus_1 * solution[i] +
1897 factor2 * matrix_diagonal_inverse[i] * (rhs[i] - tmp_vector[i]);
1904 for (std::size_t i = begin; i <
end; ++i)
1906 factor1_plus_1 * solution[i] - factor1 * solution_old[i] +
1907 factor2 * matrix_diagonal_inverse[i] * (rhs[i] - tmp_vector[i]);
1912 const Number * matrix_diagonal_inverse;
1913 const unsigned int iteration_index;
1914 const Number factor1;
1915 const Number factor2;
1917 mutable Number * tmp_vector;
1918 mutable Number * solution;
1921 template <
typename Number>
1924 VectorUpdatesRange(
const VectorUpdater<Number> &updater,
1925 const std::size_t size)
1929 VectorUpdatesRange::apply_to_subrange(0, size);
1937 ~VectorUpdatesRange()
override =
default;
1940 apply_to_subrange(
const std::size_t
begin,
1941 const std::size_t
end)
const override 1943 updater.apply_to_subrange(begin, end);
1946 const VectorUpdater<Number> &updater;
1950 template <
typename Number>
1952 vector_updates(const ::Vector<Number> & rhs,
1954 const unsigned int iteration_index,
1955 const double factor1,
1956 const double factor2,
1957 ::Vector<Number> &solution_old,
1958 ::Vector<Number> &temp_vector1,
1960 ::Vector<Number> &solution)
1962 VectorUpdater<Number> upd(rhs.begin(),
1963 jacobi.get_vector().begin(),
1967 solution_old.begin(),
1968 temp_vector1.begin(),
1970 VectorUpdatesRange<Number>(upd, rhs.size());
1973 if (iteration_index == 0)
1978 else if (iteration_index == 1)
1980 solution.swap(temp_vector1);
1981 solution_old.swap(temp_vector1);
1984 solution.swap(solution_old);
1988 template <
typename Number>
1994 const unsigned int iteration_index,
1995 const double factor1,
1996 const double factor2,
2004 VectorUpdater<Number> upd(rhs.
begin(),
2005 jacobi.get_vector().begin(),
2009 solution_old.
begin(),
2010 temp_vector1.
begin(),
2012 VectorUpdatesRange<Number>(upd, rhs.
local_size());
2015 if (iteration_index == 0)
2020 else if (iteration_index == 1)
2022 solution.
swap(temp_vector1);
2023 solution_old.
swap(temp_vector1);
2026 solution.
swap(solution_old);
2029 template <
typename MatrixType,
typename PreconditionerType>
2031 initialize_preconditioner(
2032 const MatrixType & matrix,
2033 std::shared_ptr<PreconditionerType> &preconditioner)
2036 (void)preconditioner;
2040 template <
typename MatrixType,
typename VectorType>
2042 initialize_preconditioner(
2043 const MatrixType & matrix,
2046 if (preconditioner.get() ==
nullptr || preconditioner->m() != matrix.m())
2048 if (preconditioner.get() ==
nullptr)
2052 preconditioner->m() == 0,
2054 "Preconditioner appears to be initialized but not sized correctly"));
2057 if (preconditioner->m() != matrix.m())
2059 preconditioner->get_vector().reinit(matrix.m());
2061 preconditioner->get_vector()(i) = 1. / matrix.el(i, i);
2066 template <
typename VectorType>
2070 vector = 1. / std::sqrt(static_cast<double>(vector.size()));
2071 if (vector.locally_owned_elements().is_element(0))
2075 template <
typename Number>
2077 set_initial_guess(::Vector<Number> &vector)
2083 for (
unsigned int i = 0; i < vector.size(); ++i)
2086 const Number mean_value = vector.mean_value();
2087 vector.add(-mean_value);
2090 template <
typename Number>
2105 for (
unsigned int i = 0; i < vector.
local_size(); ++i)
2108 const Number mean_value = vector.
mean_value();
2109 vector.
add(-mean_value);
2112 template <
typename Number>
2117 for (
unsigned int block = 0; block < vector.
n_blocks(); ++block)
2118 set_initial_guess(vector.
block(block));
2122 # ifdef DEAL_II_COMPILER_CUDA_AWARE 2123 template <
typename Number>
2126 const unsigned int local_size,
2130 const unsigned int index = threadIdx.x + blockDim.x * blockIdx.x;
2131 if (index < local_size)
2132 values[index] = (index + offset) % 11;
2135 template <
typename Number>
2151 const auto n_local_elements = vector.
local_size();
2154 set_initial_guess_kernel<<<n_blocks, CUDAWrappers::block_size>>>(
2155 first_local_range, n_local_elements, vector.
get_values());
2158 const Number mean_value = vector.
mean_value();
2159 vector.
add(-mean_value);
2161 # endif // DEAL_II_COMPILER_CUDA_AWARE 2163 struct EigenvalueTracker
2172 std::vector<double>
values;
2179 template <
typename MatrixType,
class VectorType,
typename PreconditionerType>
2182 const double smoothing_range,
2183 const unsigned int eig_cg_n_iterations,
2184 const double eig_cg_residual,
2185 const double max_eigenvalue)
2187 , smoothing_range(smoothing_range)
2188 , eig_cg_n_iterations(eig_cg_n_iterations)
2189 , eig_cg_residual(eig_cg_residual)
2190 , max_eigenvalue(max_eigenvalue)
2195 template <
typename MatrixType,
class VectorType,
typename PreconditionerType>
2202 degree = other_data.
degree;
2215 template <
typename MatrixType,
typename VectorType,
typename PreconditionerType>
2223 std::is_same<size_type, typename VectorType::size_type>::value,
2224 "PreconditionChebyshev and VectorType must have the same size_type.");
2229 template <
typename MatrixType,
typename VectorType,
typename PreconditionerType>
2232 const MatrixType & matrix,
2236 data = additional_data;
2238 ExcMessage(
"The degree of the Chebyshev method must be positive."));
2239 internal::PreconditionChebyshevImplementation::initialize_preconditioner(
2246 template <
typename MatrixType,
typename VectorType,
typename PreconditionerType>
2254 VectorType empty_vector;
2264 template <
typename MatrixType,
typename VectorType,
typename PreconditionerType>
2275 EigenvalueInformation info{};
2284 "Need to set at least two iterations to find eigenvalues."));
2295 internal::PreconditionChebyshevImplementation::EigenvalueTracker
2299 [&eigenvalue_tracker](
const std::vector<double> &
eigenvalues) {
2300 eigenvalue_tracker.slot(eigenvalues);
2306 internal::PreconditionChebyshevImplementation::set_initial_guess(
2321 if (eigenvalue_tracker.values.empty())
2322 info.min_eigenvalue_estimate = info.max_eigenvalue_estimate = 1.;
2325 info.min_eigenvalue_estimate = eigenvalue_tracker.values.front();
2329 info.max_eigenvalue_estimate = 1.2 * eigenvalue_tracker.values.back();
2332 info.cg_iterations = control.last_step();
2342 std::min(0.9 * info.max_eigenvalue_estimate,
2343 info.min_eigenvalue_estimate));
2352 const double actual_range = info.max_eigenvalue_estimate / alpha;
2353 const double sigma = (1. - std::sqrt(1. / actual_range)) /
2354 (1. + std::sqrt(1. / actual_range));
2360 1 + static_cast<unsigned int>(
2361 std::log(1. / eps + std::sqrt(1. / eps / eps - 1.)) /
2362 std::log(1. / sigma));
2369 ->
delta = (info.max_eigenvalue_estimate - alpha) * 0.5;
2372 ->
theta = (info.max_eigenvalue_estimate + alpha) * 0.5;
2376 using NumberType =
typename VectorType::value_type;
2379 (std::is_same<VectorType, ::Vector<NumberType>>::value ==
false &&
2380 ((std::is_same<VectorType,
2382 Vector<NumberType, MemorySpace::Host>>::value ==
2384 (std::is_same<VectorType,
2385 LinearAlgebra::distributed::
2386 Vector<NumberType, MemorySpace::CUDA>>::value ==
2391 VectorType empty_vector;
2404 template <
typename MatrixType,
typename VectorType,
typename PreconditionerType>
2410 std::lock_guard<std::mutex> lock(
mutex);
2414 internal::PreconditionChebyshevImplementation::vector_updates(
2431 for (
unsigned int k = 0; k <
data.
degree - 1; ++k)
2434 const double rhokp = 1. / (2. * sigma - rhok);
2435 const double factor1 = rhokp * rhok, factor2 = 2. * rhokp /
delta;
2437 internal::PreconditionChebyshevImplementation::vector_updates(
2452 template <
typename MatrixType,
typename VectorType,
typename PreconditionerType>
2458 std::lock_guard<std::mutex> lock(
mutex);
2462 internal::PreconditionChebyshevImplementation::vector_updates(
2477 for (
unsigned int k = 0; k <
data.
degree - 1; ++k)
2480 const double rhokp = 1. / (2. * sigma - rhok);
2481 const double factor1 = rhokp * rhok, factor2 = 2. * rhokp /
delta;
2483 internal::PreconditionChebyshevImplementation::vector_updates(
2498 template <
typename MatrixType,
typename VectorType,
typename PreconditionerType>
2504 std::lock_guard<std::mutex> lock(
mutex);
2509 internal::PreconditionChebyshevImplementation::vector_updates(
2524 for (
unsigned int k = 0; k <
data.
degree - 1; ++k)
2527 const double rhokp = 1. / (2. * sigma - rhok);
2528 const double factor1 = rhokp * rhok, factor2 = 2. * rhokp /
delta;
2530 internal::PreconditionChebyshevImplementation::vector_updates(
2545 template <
typename MatrixType,
typename VectorType,
typename PreconditionerType>
2551 std::lock_guard<std::mutex> lock(
mutex);
2556 internal::PreconditionChebyshevImplementation::vector_updates(
2571 for (
unsigned int k = 0; k <
data.
degree - 1; ++k)
2574 const double rhokp = 1. / (2. * sigma - rhok);
2575 const double factor1 = rhokp * rhok, factor2 = 2. * rhokp /
delta;
2577 internal::PreconditionChebyshevImplementation::vector_updates(
2592 template <
typename MatrixType,
typename VectorType,
typename PreconditionerType>
2604 template <
typename MatrixType,
typename VectorType,
typename PreconditionerType>
void Tvmult(VectorType &, const VectorType &) const
size_type local_size() const
static const unsigned int invalid_unsigned_int
void set_zero(VectorType &vec) const
void initialize(const MatrixType &A, const typename BaseClass::AdditionalData ¶meters=typename BaseClass::AdditionalData())
const_iterator end() const
Contents is actually a matrix.
types::global_dof_index size_type
void Tstep(VectorType &x, const VectorType &rhs) const
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
void vmult(VectorType &dst, const VectorType &src) const
void Tvmult_add(VectorType &, const VectorType &) const
PreconditionRelaxation< MatrixType >::AdditionalData parameters
void vmult(VectorType &, const VectorType &) const
Number * get_values() const
const MatrixType & matrix
AdditionalData(const double relaxation=1.)
void(MatrixType::*)(VectorType &, const VectorType &) const function_ptr
static ::ExceptionBase & ExcNotInitialized()
void initialize(const MatrixType &matrix, const AdditionalData &additional_data=AdditionalData())
std::shared_ptr< PreconditionerType > preconditioner
#define AssertThrow(cond, exc)
AdditionalData(const unsigned int degree=1, const double smoothing_range=0., const unsigned int eig_cg_n_iterations=8, const double eig_cg_residual=1e-2, const double max_eigenvalue=1)
AdditionalData(const double relaxation=1.)
#define AssertCudaKernel()
const std::vector< size_type > * permutation
SymmetricTensor< 2, dim, Number > epsilon(const Tensor< 2, dim, Number > &Grad_u)
const function_ptr precondition
void vmult(VectorType &, const VectorType &) const
void initialize(const AdditionalData ¶meters)
Number local_element(const size_type local_index) const
void Tvmult(VectorType &, const VectorType &) const
static ::ExceptionBase & ExcMessage(std::string arg1)
typename MatrixType::size_type size_type
AdditionalData & operator=(const AdditionalData &other_data)
void Tvmult(VectorType &dst, const VectorType &src) const
void vmult(VectorType &, const VectorType &) const
const std::vector< size_type > & inverse_permutation
#define Assert(cond, exc)
typename PreconditionRelaxation< MatrixType >::AdditionalData AdditionalData
bool eigenvalues_are_initialized
void step(VectorType &x, const VectorType &rhs) const
#define DEAL_II_NAMESPACE_CLOSE
VectorType::value_type * end(VectorType &V)
typename MatrixType::size_type size_type
void vmult_add(VectorType &, const VectorType &) const
std::enable_if< IsBlockVector< VectorType >::value, unsigned int >::type n_blocks(const VectorType &vector)
void Tvmult(VectorType &, const VectorType &) const
void initialize(const MatrixType &A, const AdditionalData ¶meters=AdditionalData())
void Tvmult(VectorType &, const VectorType &) const
void Tstep(VectorType &x, const VectorType &rhs) const
EigenvalueInformation estimate_eigenvalues(const VectorType &src) const
void vmult(VectorType &, const VectorType &) const
void Tstep(VectorType &x, const VectorType &rhs) const
void solve(const MatrixType &A, VectorType &x, const VectorType &b, const PreconditionerType &preconditioner)
void vmult(VectorType &dst, const VectorType &src) const
AdditionalData(const std::vector< size_type > &permutation, const std::vector< size_type > &inverse_permutation, const typename PreconditionRelaxation< MatrixType >::AdditionalData ¶meters=typename PreconditionRelaxation< MatrixType >::AdditionalData())
void Tvmult_add(VectorType &, const VectorType &) const
void Tvmult(VectorType &, const VectorType &) const
void vmult(VectorType &, const VectorType &) const
unsigned int global_dof_index
virtual ::IndexSet locally_owned_elements() const override
SmartPointer< const MatrixType, PreconditionRelaxation< MatrixType > > A
void Tvmult(VectorType &, const VectorType &) const
void step(VectorType &x, const VectorType &rhs) const
AffineConstraints< double > constraints
void initialize(const MatrixType &matrix, const AdditionalData &additional_data=AdditionalData())
boost::signals2::connection connect_eigenvalues_slot(const std::function< void(const std::vector< double > &)> &slot, const bool every_iteration=false)
unsigned int n_blocks() const
void swap(Vector< Number, MemorySpace > &v)
ARKode< VectorType > * solver
#define DEAL_II_NAMESPACE_OPEN
VectorType::value_type * begin(VectorType &V)
T min(const T &t, const MPI_Comm &mpi_communicator)
typename MatrixType::size_type size_type
const std::vector< size_type > * inverse_permutation
void Tstep(VectorType &dst, const VectorType &src) const
void step(VectorType &dst, const VectorType &src) const
void step(VectorType &x, const VectorType &rhs) const
unsigned int minimum_parallel_grain_size
virtual void add(const Number a) override
unsigned int eig_cg_n_iterations
const_iterator begin() const
typename PreconditionRelaxation< MatrixType >::AdditionalData AdditionalData
void vmult(VectorType &, const VectorType &) const
const std::vector< size_type > & permutation
PreconditionUseMatrix(const MatrixType &M, const function_ptr method)
typename PreconditionRelaxation< MatrixType >::AdditionalData AdditionalData
Eigenvalue vector is filled.
virtual Number mean_value() const override
void vmult_add(VectorType &, const VectorType &) const
BlockType & block(const unsigned int i)
std::vector< std::size_t > pos_right_of_diagonal
size_type nth_index_in_set(const size_type local_index) const
T max(const T &t, const MPI_Comm &mpi_communicator)
#define DEAL_II_OPENMP_SIMD_PRAGMA
SmartPointer< const MatrixType, PreconditionChebyshev< MatrixType, VectorType, PreconditionerType > > matrix_ptr
static ::ExceptionBase & ExcInternalError()
void initialize(const MatrixType &A, const std::vector< size_type > &permutation, const std::vector< size_type > &inverse_permutation, const typename PreconditionRelaxation< MatrixType >::AdditionalData ¶meters=typename PreconditionRelaxation< MatrixType >::AdditionalData())