15#ifndef dealii_solver_cg_h
16#define dealii_solver_cg_h
31#include <boost/signals2.hpp>
44 template <
typename,
typename>
196template <
typename VectorType = Vector<
double>>
256 template <
typename MatrixType,
typename PreconditionerType>
260 void solve(const MatrixType &A,
263 const PreconditionerType &preconditioner);
271 boost::signals2::connection
272 connect_coefficients_slot(
273 const
std::function<
void(typename VectorType::value_type,
274 typename VectorType::value_type)> &slot);
282 boost::signals2::connection
283 connect_condition_number_slot(const
std::function<
void(
double)> &slot,
284 const
bool every_iteration = false);
292 boost::signals2::connection
293 connect_eigenvalues_slot(
294 const
std::function<
void(const
std::vector<
double> &)> &slot,
295 const
bool every_iteration = false);
304 print_vectors(const
unsigned int step,
307 const VectorType &d) const;
315 compute_eigs_and_cond(
316 const
std::vector<typename VectorType::value_type> &diagonal,
317 const
std::vector<typename VectorType::value_type> &offdiagonal,
318 const
boost::signals2::signal<
void(const
std::vector<
double> &)>
320 const
boost::signals2::signal<
void(
double)> &cond_signal);
330 boost::signals2::signal<
void(typename VectorType::value_type,
331 typename VectorType::value_type)>
338 boost::signals2::signal<
void(
double)> condition_number_signal;
344 boost::signals2::signal<
void(
double)> all_condition_numbers_signal;
350 boost::signals2::signal<
void(const
std::vector<
double> &)> eigenvalues_signal;
356 boost::signals2::signal<
void(const
std::vector<
double> &)>
357 all_eigenvalues_signal;
367 bool determine_beta_by_flexible_formula;
397template <typename VectorType =
Vector<
double>>
465template <
typename VectorType>
468 const bool use_default_residual)
469 : use_default_residual(use_default_residual)
472template <
typename VectorType>
476 const AdditionalData &
data)
478 , additional_data(
data)
479 , determine_beta_by_flexible_formula(false)
484template <
typename VectorType>
488 , additional_data(
data)
489 , determine_beta_by_flexible_formula(false)
494template <
typename VectorType>
499 const VectorType &)
const
504template <
typename VectorType>
507 const std::vector<typename VectorType::value_type> &diagonal,
508 const std::vector<typename VectorType::value_type> &offdiagonal,
509 const boost::signals2::signal<
void(
const std::vector<double> &)>
511 const boost::signals2::signal<
void(
double)> &cond_signal)
514 if (!cond_signal.empty() || !eigenvalues_signal.empty())
518 for (size_type i = 0; i <
diagonal.size(); ++i)
522 T(i, i + 1) = offdiagonal[i];
524 T.compute_eigenvalues();
528 auto condition_number =
T.eigenvalue(
T.n() - 1) /
T.eigenvalue(0);
531 cond_signal(
std::abs(condition_number));
535 if (!eigenvalues_signal.empty())
538 for (
unsigned int j = 0; j <
T.n(); ++j)
564 typename PreconditionerType>
565 struct IterationWorkerBase
567 using Number =
typename VectorType::value_type;
570 const PreconditionerType &preconditioner;
593 Number r_dot_preconditioner_dot_r;
596 double residual_norm;
597 Number previous_alpha;
598 const bool use_default_residual;
601 IterationWorkerBase(
const MatrixType &A,
602 const PreconditionerType &preconditioner,
607 const bool &use_default_residual)
609 , preconditioner(preconditioner)
617 , explicit_r_pointer(memory)
622 , explicit_r(*explicit_r_pointer)
623 , r_dot_preconditioner_dot_r(Number())
627 , previous_alpha(Number())
628 , use_default_residual(use_default_residual)
641 if (!use_default_residual)
642 explicit_r.reinit(x,
true);
654 residual_norm = r.l2_norm();
664 typename PreconditionerType,
666 struct IterationWorker
667 :
public IterationWorkerBase<VectorType, MatrixType, PreconditionerType>
670 IterationWorkerBase<VectorType, MatrixType, PreconditionerType>;
673 IterationWorker(
const MatrixType &A,
674 const PreconditionerType &preconditioner,
679 const bool &use_default_residual)
686 use_default_residual)
690 using BaseClass::alpha;
692 using BaseClass::beta;
693 using BaseClass::explicit_r;
695 using BaseClass::preconditioner;
697 using BaseClass::r_dot_preconditioner_dot_r;
698 using BaseClass::residual_norm;
699 using BaseClass::use_default_residual;
705 do_iteration(
const unsigned int iteration_index)
707 using Number =
typename VectorType::value_type;
709 const Number previous_r_dot_preconditioner_dot_r =
710 r_dot_preconditioner_dot_r;
712 if (std::is_same_v<PreconditionerType, PreconditionIdentity> ==
false)
714 preconditioner.vmult(v, r);
715 r_dot_preconditioner_dot_r = r * v;
718 r_dot_preconditioner_dot_r = residual_norm * residual_norm;
721 std::is_same_v<PreconditionerType, PreconditionIdentity> ? r : v;
723 if (iteration_index > 1)
728 r_dot_preconditioner_dot_r / previous_r_dot_preconditioner_dot_r;
730 beta -= (r * z) / previous_r_dot_preconditioner_dot_r;
731 p.sadd(beta, 1., direction);
734 p.equ(1., direction);
741 const Number p_dot_A_dot_p = p * v;
744 this->previous_alpha = alpha;
745 alpha = r_dot_preconditioner_dot_r / p_dot_A_dot_p;
750 if (use_default_residual)
761 A.vmult(explicit_r, x);
762 explicit_r.add(-1, b);
763 residual_norm = explicit_r.l2_norm();
768 finalize_after_convergence(
const unsigned int)
780 template <
typename MatrixType,
typename VectorType>
781 using vmult_functions_t =
decltype(std::declval<const MatrixType>().vmult(
782 std::declval<VectorType &>(),
783 std::declval<const VectorType &>(),
785 const std::function<
void(
const unsigned int,
const unsigned int)> &>(),
786 std::declval<
const std::function<
void(
const unsigned int,
787 const unsigned int)> &>()));
789 template <
typename MatrixType,
typename VectorType>
790 constexpr bool has_vmult_functions =
791 is_supported_operation<vmult_functions_t, MatrixType, VectorType>;
796 template <
typename PreconditionerType>
797 using apply_to_subrange_t =
798 decltype(std::declval<const PreconditionerType>()
799 .apply_to_subrange(0U, 0U,
nullptr,
nullptr));
801 template <
typename PreconditionerType>
802 constexpr bool has_apply_to_subrange =
803 is_supported_operation<apply_to_subrange_t, PreconditionerType>;
808 template <
typename PreconditionerType>
810 decltype(std::declval<const PreconditionerType>().apply(0U, 0.0));
812 template <
typename PreconditionerType>
813 constexpr bool has_apply =
814 is_supported_operation<apply_t, PreconditionerType>;
822 typename PreconditionerType>
823 struct IterationWorker<
827 std::enable_if_t<has_vmult_functions<MatrixType, VectorType> &&
828 (has_apply_to_subrange<PreconditionerType> ||
829 has_apply<PreconditionerType>)&&std::
830 is_same_v<VectorType,
831 LinearAlgebra::distributed::Vector<
832 typename VectorType::value_type,
835 :
public IterationWorkerBase<VectorType, MatrixType, PreconditionerType>
837 using Number =
typename VectorType::value_type;
839 Number next_r_dot_preconditioner_dot_r;
840 Number previous_beta;
842 IterationWorker(
const MatrixType &A,
843 const PreconditionerType &preconditioner,
848 const bool &use_default_residual)
856 use_default_residual)
857 , next_r_dot_preconditioner_dot_r(0.)
864 do_iteration(
const unsigned int iteration_index)
866 if (iteration_index > 1)
868 previous_beta = this->beta;
869 this->beta = next_r_dot_preconditioner_dot_r /
870 this->r_dot_preconditioner_dot_r;
873 std::array<VectorizedArray<Number>, 7> vectorized_sums = {};
878 [&](
const unsigned int begin,
const unsigned int end) {
879 operation_before_loop(iteration_index, begin, end);
881 [&](
const unsigned int begin,
const unsigned int end) {
882 operation_after_loop(begin, end, vectorized_sums);
885 std::array<Number, 7> scalar_sums;
886 for (
unsigned int i = 0; i < 7; ++i)
887 scalar_sums[i] = vectorized_sums[i][0];
888 for (
unsigned int l = 1; l < VectorizedArray<Number>::size(); ++
l)
889 for (
unsigned int i = 0; i < 7; ++i)
890 scalar_sums[i] += vectorized_sums[i][l];
894 this->r.get_mpi_communicator(),
897 this->r_dot_preconditioner_dot_r = scalar_sums[6];
899 const Number p_dot_A_dot_p = scalar_sums[0];
902 this->previous_alpha = this->alpha;
903 this->alpha = this->r_dot_preconditioner_dot_r / p_dot_A_dot_p;
909 this->alpha * (-2. * scalar_sums[2] + this->alpha * scalar_sums[1])));
911 next_r_dot_preconditioner_dot_r =
std::abs(
912 this->r_dot_preconditioner_dot_r +
913 this->alpha * (-2. * scalar_sums[4] + this->alpha * scalar_sums[5]));
918 template <
typename U =
void>
919 std::enable_if_t<has_apply<PreconditionerType>,
U>
920 operation_before_loop(
const unsigned int iteration_index,
921 const unsigned int start_range,
922 const unsigned int end_range)
const
924 Number *x = this->x.begin();
925 Number *r = this->r.begin();
926 Number *p = this->p.begin();
927 Number *v = this->v.begin();
928 const Number alpha = this->alpha;
929 const Number beta = this->beta;
931 const unsigned int end_regular =
932 start_range + (end_range - start_range) / n_lanes * n_lanes;
933 if (iteration_index == 1)
938 for (
unsigned int j = start_range; j < end_regular; j += n_lanes)
943 for (
unsigned int l = 0;
l < n_lanes; ++
l)
944 pj[l] = this->preconditioner.apply(j + l, rj[l]);
949 for (
unsigned int j = end_regular; j < end_range; ++j)
951 p[j] = this->preconditioner.apply(j, r[j]);
955 else if (iteration_index % 2 == 0 && beta != Number())
957 for (
unsigned int j = start_range; j < end_regular; j += n_lanes)
965 for (
unsigned int l = 0;
l < n_lanes; ++
l)
966 prec_rj[l] = this->preconditioner.apply(j + l, rj[l]);
968 pj = beta * pj + prec_rj;
973 for (
unsigned int j = end_regular; j < end_range; ++j)
975 r[j] -= alpha * v[j];
976 p[j] = beta * p[j] + this->preconditioner.apply(j, r[j]);
980 else if (iteration_index % 2 == 0 && beta == Number())
985 for (
unsigned int j = start_range; j < end_regular; j += n_lanes)
997 for (
unsigned int l = 0;
l < n_lanes; ++
l)
998 prec_rj[l] = this->preconditioner.apply(j + l, rj[l]);
999 prec_rj.
store(p + j);
1003 for (
unsigned int j = end_regular; j < end_range; ++j)
1005 r[j] -= alpha * v[j];
1006 x[j] += alpha * p[j];
1007 p[j] = this->preconditioner.apply(j, r[j]);
1013 const Number alpha_plus_previous_alpha_over_beta =
1014 alpha + this->previous_alpha / this->previous_beta;
1015 const Number previous_alpha_over_beta =
1016 this->previous_alpha / this->previous_beta;
1017 for (
unsigned int j = start_range; j < end_regular; j += n_lanes)
1022 xj += alpha_plus_previous_alpha_over_beta * pj;
1026 for (
unsigned int l = 0;
l < n_lanes; ++
l)
1028 prec_rj[
l] = this->preconditioner.apply(j + l, rj[l]);
1029 prec_vj[
l] = this->preconditioner.apply(j + l, vj[l]);
1031 xj -= previous_alpha_over_beta * prec_rj;
1035 prec_rj -= alpha * prec_vj;
1036 pj = beta * pj + prec_rj;
1041 for (
unsigned int j = end_regular; j < end_range; ++j)
1043 x[j] += alpha_plus_previous_alpha_over_beta * p[j];
1044 x[j] -= previous_alpha_over_beta *
1045 this->preconditioner.apply(j, r[j]);
1046 r[j] -= alpha * v[j];
1047 p[j] = beta * p[j] + this->preconditioner.apply(j, r[j]);
1055 template <
typename U =
void>
1056 std::enable_if_t<has_apply<PreconditionerType>,
U>
1057 operation_after_loop(
1058 const unsigned int start_range,
1059 const unsigned int end_range,
1062 const Number *r = this->r.begin();
1063 const Number *p = this->p.begin();
1064 const Number *v = this->v.begin();
1065 std::array<VectorizedArray<Number>, 7> my_sums = {};
1067 const unsigned int end_regular =
1068 start_range + (end_range - start_range) / n_lanes * n_lanes;
1069 for (
unsigned int j = start_range; j < end_regular; j += n_lanes)
1076 for (
unsigned int l = 0;
l < n_lanes; ++
l)
1078 prec_vj[
l] = this->preconditioner.apply(j + l, vj[l]);
1079 prec_rj[
l] = this->preconditioner.apply(j + l, rj[l]);
1081 my_sums[0] += pj * vj;
1082 my_sums[1] += vj * vj;
1083 my_sums[2] += rj * vj;
1084 my_sums[3] += rj * rj;
1085 my_sums[4] += rj * prec_vj;
1086 my_sums[5] += vj * prec_vj;
1087 my_sums[6] += rj * prec_rj;
1089 for (
unsigned int j = end_regular; j < end_range; ++j)
1091 const Number prec_v = this->preconditioner.apply(j, v[j]);
1092 const Number prec_r = this->preconditioner.apply(j, r[j]);
1093 my_sums[0][0] += p[j] * v[j];
1094 my_sums[1][0] += v[j] * v[j];
1095 my_sums[2][0] += r[j] * v[j];
1096 my_sums[3][0] += r[j] * r[j];
1097 my_sums[4][0] += r[j] * prec_v;
1098 my_sums[5][0] += v[j] * prec_v;
1099 my_sums[6][0] += r[j] * prec_r;
1101 for (
unsigned int i = 0; i < vectorized_sums.size(); ++i)
1102 vectorized_sums[i] += my_sums[i];
1107 template <
typename U =
void>
1108 std::enable_if_t<has_apply<PreconditionerType>,
U>
1109 finalize_after_convergence(
const unsigned int iteration_index)
1111 if (iteration_index % 2 == 1 || this->beta == Number())
1112 this->x.add(this->alpha, this->p);
1115 using Number =
typename VectorType::value_type;
1116 const unsigned int end_range = this->x.locally_owned_size();
1118 Number *x = this->x.begin();
1119 Number *r = this->r.begin();
1120 Number *p = this->p.begin();
1126 const Number alpha_plus_previous_alpha_over_beta =
1127 this->alpha + this->previous_alpha / this->beta;
1128 const Number previous_alpha_over_beta =
1129 this->previous_alpha / this->beta;
1132 for (
unsigned int j = 0; j < end_range; ++j)
1134 x[j] += alpha_plus_previous_alpha_over_beta * p[j] -
1135 previous_alpha_over_beta *
1136 this->preconditioner.apply(j, r[j]);
1144 template <
typename U =
void>
1145 std::enable_if_t<!has_apply<PreconditionerType>,
U>
1146 operation_before_loop(
const unsigned int iteration_index,
1147 const unsigned int start_range,
1148 const unsigned int end_range)
const
1150 Number *x = this->x.begin() + start_range;
1151 Number *r = this->r.begin() + start_range;
1152 Number *p = this->p.begin() + start_range;
1153 Number *v = this->v.begin() + start_range;
1154 const Number alpha = this->alpha;
1155 const Number beta = this->beta;
1156 constexpr unsigned int grain_size = 128;
1157 std::array<Number, grain_size> prec_r;
1158 if (iteration_index == 1)
1160 for (
unsigned int j = start_range; j < end_range; j += grain_size)
1162 const unsigned int length =
std::min(grain_size, end_range - j);
1163 this->preconditioner.apply_to_subrange(j,
1168 for (
unsigned int i = 0; i < length; ++i)
1178 else if (iteration_index % 2 == 0 && beta != Number())
1180 for (
unsigned int j = start_range; j < end_range; j += grain_size)
1182 const unsigned int length =
std::min(grain_size, end_range - j);
1184 for (
unsigned int i = 0; i < length; ++i)
1185 r[i] -= this->alpha * v[i];
1186 this->preconditioner.apply_to_subrange(j,
1191 for (
unsigned int i = 0; i < length; ++i)
1193 p[i] = this->beta * p[i] + prec_r[i];
1201 else if (iteration_index % 2 == 0 && beta == Number())
1206 for (
unsigned int j = start_range; j < end_range; j += grain_size)
1208 const unsigned int length =
std::min(grain_size, end_range - j);
1210 for (
unsigned int i = 0; i < length; ++i)
1211 r[i] -= this->alpha * v[i];
1212 this->preconditioner.apply_to_subrange(j,
1217 for (
unsigned int i = 0; i < length; ++i)
1219 x[i] += this->alpha * p[i];
1230 const Number alpha_plus_previous_alpha_over_beta =
1231 this->alpha + this->previous_alpha / this->previous_beta;
1232 const Number previous_alpha_over_beta =
1233 this->previous_alpha / this->previous_beta;
1234 for (
unsigned int j = start_range; j < end_range; j += grain_size)
1236 const unsigned int length =
std::min(grain_size, end_range - j);
1237 this->preconditioner.apply_to_subrange(j,
1242 for (
unsigned int i = 0; i < length; ++i)
1244 x[i] += alpha_plus_previous_alpha_over_beta * p[i] -
1245 previous_alpha_over_beta * prec_r[i];
1246 r[i] -= this->alpha * v[i];
1248 this->preconditioner.apply_to_subrange(j,
1253 for (
unsigned int i = 0; i < length; ++i)
1255 p[i] = this->beta * p[i] + prec_r[i];
1269 template <
typename U =
void>
1270 std::enable_if_t<!has_apply<PreconditionerType>,
U>
1271 operation_after_loop(
1272 const unsigned int start_range,
1273 const unsigned int end_range,
1276 const Number *r = this->r.begin();
1277 const Number *p = this->p.begin();
1278 const Number *v = this->v.begin();
1279 std::array<VectorizedArray<Number>, 7> my_sums = {};
1280 constexpr unsigned int grain_size = 128;
1283 const unsigned int end_regular =
1284 start_range + (end_range - start_range) / grain_size * grain_size;
1285 std::array<Number, grain_size> prec_r;
1286 std::array<Number, grain_size> prec_v;
1287 for (
unsigned int j = start_range; j < end_regular; j += grain_size)
1289 this->preconditioner.apply_to_subrange(j,
1293 this->preconditioner.apply_to_subrange(j,
1298 for (
unsigned int i = 0; i < grain_size;
1304 prec_rj.
load(prec_r.data() + i);
1305 prec_vj.
load(prec_v.data() + i);
1307 my_sums[0] += pj * vj;
1308 my_sums[1] += vj * vj;
1309 my_sums[2] += rj * vj;
1310 my_sums[3] += rj * rj;
1311 my_sums[4] += rj * prec_vj;
1312 my_sums[5] += vj * prec_vj;
1313 my_sums[6] += rj * prec_rj;
1316 const unsigned int length = end_range - end_regular;
1318 this->preconditioner.apply_to_subrange(end_regular,
1319 end_regular + length,
1322 this->preconditioner.apply_to_subrange(end_regular,
1323 end_regular + length,
1326 for (
unsigned int j = end_regular; j < end_range; ++j)
1328 my_sums[0][0] += p[j] * v[j];
1329 my_sums[1][0] += v[j] * v[j];
1330 my_sums[2][0] += r[j] * v[j];
1331 my_sums[3][0] += r[j] * r[j];
1332 my_sums[4][0] += r[j] * prec_v[j - end_regular];
1333 my_sums[5][0] += v[j] * prec_v[j - end_regular];
1334 my_sums[6][0] += r[j] * prec_r[j - end_regular];
1336 for (
unsigned int i = 0; i < vectorized_sums.size(); ++i)
1337 vectorized_sums[i] += my_sums[i];
1343 template <
typename U =
void>
1344 std::enable_if_t<!has_apply<PreconditionerType>,
U>
1345 finalize_after_convergence(
const unsigned int iteration_index)
1347 if (iteration_index % 2 == 1 || this->beta == Number())
1348 this->x.add(this->alpha, this->p);
1351 const unsigned int end_range = this->x.locally_owned_size();
1353 Number *x = this->x.begin();
1354 Number *r = this->r.begin();
1355 Number *p = this->p.begin();
1356 const Number alpha_plus_previous_alpha_over_beta =
1357 this->alpha + this->previous_alpha / this->beta;
1358 const Number previous_alpha_over_beta =
1359 this->previous_alpha / this->beta;
1361 constexpr unsigned int grain_size = 128;
1362 std::array<Number, grain_size> prec_r;
1363 for (
unsigned int j = 0; j < end_range; j += grain_size)
1365 const unsigned int length =
std::min(grain_size, end_range - j);
1366 this->preconditioner.apply_to_subrange(j,
1371 for (
unsigned int i = 0; i < length; ++i)
1372 x[i] += alpha_plus_previous_alpha_over_beta * p[i] -
1373 previous_alpha_over_beta * prec_r[i];
1387template <
typename VectorType>
1389template <
typename MatrixType,
typename PreconditionerType>
1395 const VectorType &b,
1396 const PreconditionerType &preconditioner)
1398 using number =
typename VectorType::value_type;
1405 const bool do_eigenvalues =
1406 !condition_number_signal.empty() || !all_condition_numbers_signal.empty() ||
1407 !eigenvalues_signal.empty() || !all_eigenvalues_signal.empty();
1410 std::vector<typename VectorType::value_type>
diagonal;
1411 std::vector<typename VectorType::value_type> offdiagonal;
1413 typename VectorType::value_type eigen_beta_alpha = 0;
1417 internal::SolverCG::
1418 IterationWorker<VectorType, MatrixType, PreconditionerType>
1421 determine_beta_by_flexible_formula,
1429 solver_state = this->iteration_status(0, worker.residual_norm, x);
1437 worker.do_iteration(it);
1439 print_vectors(it, x, worker.r, worker.p);
1443 this->coefficients_signal(worker.previous_alpha, worker.beta);
1448 diagonal.push_back(number(1.) / worker.previous_alpha +
1450 eigen_beta_alpha = worker.beta / worker.previous_alpha;
1451 offdiagonal.push_back(
std::sqrt(worker.beta) /
1452 worker.previous_alpha);
1454 compute_eigs_and_cond(diagonal,
1456 all_eigenvalues_signal,
1457 all_condition_numbers_signal);
1460 solver_state = this->iteration_status(it, worker.residual_norm, x);
1463 worker.finalize_after_convergence(it);
1465 compute_eigs_and_cond(diagonal,
1468 condition_number_signal);
1476template <
typename VectorType>
1479 const std::function<
void(
typename VectorType::value_type,
1480 typename VectorType::value_type)> &slot)
1482 return coefficients_signal.connect(slot);
1487template <
typename VectorType>
1490 const std::function<
void(
double)> &slot,
1491 const bool every_iteration)
1493 if (every_iteration)
1495 return all_condition_numbers_signal.connect(slot);
1499 return condition_number_signal.connect(slot);
1505template <
typename VectorType>
1508 const std::function<
void(
const std::vector<double> &)> &slot,
1509 const bool every_iteration)
1511 if (every_iteration)
1513 return all_eigenvalues_signal.connect(slot);
1517 return eigenvalues_signal.connect(slot);
1523template <
typename VectorType>
1526 const bool use_default_residual)
1527 : use_default_residual(use_default_residual)
1531template <
typename VectorType>
1535 const AdditionalData &
data)
1538 this->determine_beta_by_flexible_formula =
true;
1539 this->additional_data =
data;
1544template <
typename VectorType>
1547 const AdditionalData &
data)
1550 this->determine_beta_by_flexible_formula =
true;
1551 this->additional_data =
data;
SolverCG(SolverControl &cn, VectorMemory< VectorType > &mem, const AdditionalData &data=AdditionalData())
virtual void print_vectors(const unsigned int step, const VectorType &x, const VectorType &r, const VectorType &d) const
boost::signals2::connection connect_condition_number_slot(const std::function< void(double)> &slot, const bool every_iteration=false)
boost::signals2::connection connect_eigenvalues_slot(const std::function< void(const std::vector< double > &)> &slot, const bool every_iteration=false)
boost::signals2::connection connect_coefficients_slot(const std::function< void(typename VectorType::value_type, typename VectorType::value_type)> &slot)
virtual ~SolverCG() override=default
void solve(const MatrixType &A, VectorType &x, const VectorType &b, const PreconditionerType &preconditioner)
SolverCG(SolverControl &cn, const AdditionalData &data=AdditionalData())
static void compute_eigs_and_cond(const std::vector< typename VectorType::value_type > &diagonal, const std::vector< typename VectorType::value_type > &offdiagonal, const boost::signals2::signal< void(const std::vector< double > &)> &eigenvalues_signal, const boost::signals2::signal< void(double)> &cond_signal)
@ iterate
Continue iteration.
@ success
Stop iteration, goal reached.
AdditionalData additional_data
SolverFlexibleCG(SolverControl &cn, const AdditionalData &data=AdditionalData())
SolverFlexibleCG(SolverControl &cn, VectorMemory< VectorType > &mem, const AdditionalData &data=AdditionalData())
static constexpr std::size_t size()
void store(OtherNumber *ptr) const
void load(const OtherNumber *ptr)
#define DEAL_II_OPENMP_SIMD_PRAGMA
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_CXX20_REQUIRES(condition)
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcDivideByZero()
#define AssertThrow(cond, exc)
std::vector< index_type > data
@ diagonal
Matrix is diagonal.
Tpetra::Vector< Number, LO, GO, NodeType< MemorySpace > > VectorType
Tpetra::CrsMatrix< Number, LO, GO, NodeType< MemorySpace > > MatrixType
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
VectorType::value_type * end(VectorType &V)
VectorType::value_type * begin(VectorType &V)
T sum(const T &t, const MPI_Comm mpi_communicator)
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > sqrt(const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > abs(const ::VectorizedArray< Number, width > &)
unsigned int global_dof_index
AdditionalData(const bool use_default_residual=true)
bool use_default_residual
bool use_default_residual
AdditionalData(const bool use_default_residual=true)
std::array< Number, 1 > eigenvalues(const SymmetricTensor< 2, 1, Number > &T)