16 #ifndef dealii_solver_cg_h
17 #define dealii_solver_cg_h
42 template <
typename,
typename>
176 template <
typename VectorType = Vector<
double>>
214 template <
typename MatrixType,
typename PreconditionerType>
219 const PreconditionerType &preconditioner);
227 boost::signals2::connection
229 const std::function<
void(
typename VectorType::value_type,
230 typename VectorType::value_type)> &slot);
238 boost::signals2::connection
240 const bool every_iteration =
false);
248 boost::signals2::connection
250 const std::function<
void(
const std::vector<double> &)> &slot,
251 const bool every_iteration =
false);
263 const VectorType &
d)
const;
272 const std::vector<typename VectorType::value_type> &
diagonal,
273 const std::vector<typename VectorType::value_type> &offdiagonal,
274 const boost::signals2::signal<
void(
const std::vector<double> &)>
276 const boost::signals2::signal<
void(
double)> &cond_signal);
286 boost::signals2::signal<void(
typename VectorType::value_type,
287 typename VectorType::value_type)>
312 boost::signals2::signal<void(
const std::vector<double> &)>
353 template <
typename VectorType = Vector<
double>>
394 template <
typename VectorType>
397 const AdditionalData &data)
399 , additional_data(data)
400 , determine_beta_by_flexible_formula(false)
405 template <
typename VectorType>
408 , additional_data(data)
409 , determine_beta_by_flexible_formula(false)
414 template <
typename VectorType>
419 const VectorType &)
const
424 template <
typename VectorType>
427 const std::vector<typename VectorType::value_type> &
diagonal,
428 const std::vector<typename VectorType::value_type> &offdiagonal,
429 const boost::signals2::signal<
void(
const std::vector<double> &)>
431 const boost::signals2::signal<
void(
double)> &cond_signal)
434 if (!cond_signal.empty() || !eigenvalues_signal.empty())
442 T(i, i + 1) = offdiagonal[i];
444 T.compute_eigenvalues();
448 auto condition_number =
T.eigenvalue(
T.n() - 1) /
T.eigenvalue(0);
451 cond_signal(std::abs(condition_number));
455 if (!eigenvalues_signal.empty())
458 for (
unsigned int j = 0; j <
T.n(); ++j)
481 template <
typename VectorType,
483 typename PreconditionerType>
484 struct IterationWorkerBase
486 using Number =
typename VectorType::value_type;
489 const PreconditionerType &preconditioner;
509 Number r_dot_preconditioner_dot_r;
512 double residual_norm;
513 Number previous_alpha;
515 IterationWorkerBase(
const MatrixType &
A,
516 const PreconditionerType &preconditioner,
521 , preconditioner(preconditioner)
532 , r_dot_preconditioner_dot_r(Number())
536 , previous_alpha(Number())
540 startup(
const VectorType &
b)
560 residual_norm = r.l2_norm();
568 template <
typename VectorType,
570 typename PreconditionerType,
572 struct IterationWorker
573 :
public IterationWorkerBase<VectorType, MatrixType, PreconditionerType>
576 IterationWorkerBase<VectorType, MatrixType, PreconditionerType>;
578 IterationWorker(
const MatrixType &
A,
579 const PreconditionerType &preconditioner,
583 : BaseClass(
A, preconditioner, flexible, memory, x)
587 using BaseClass::alpha;
588 using BaseClass::beta;
590 using BaseClass::preconditioner;
592 using BaseClass::r_dot_preconditioner_dot_r;
593 using BaseClass::residual_norm;
599 do_iteration(
const unsigned int iteration_index)
601 using Number =
typename VectorType::value_type;
603 const Number previous_r_dot_preconditioner_dot_r =
604 r_dot_preconditioner_dot_r;
606 if (std::is_same_v<PreconditionerType, PreconditionIdentity> ==
false)
608 preconditioner.vmult(v, r);
609 r_dot_preconditioner_dot_r = r * v;
612 r_dot_preconditioner_dot_r = residual_norm * residual_norm;
614 const VectorType &direction =
615 std::is_same_v<PreconditionerType, PreconditionIdentity> ? r : v;
617 if (iteration_index > 1)
619 Assert(std::abs(previous_r_dot_preconditioner_dot_r) != 0.,
622 r_dot_preconditioner_dot_r / previous_r_dot_preconditioner_dot_r;
624 beta -= (r * z) / previous_r_dot_preconditioner_dot_r;
625 p.sadd(beta, 1., direction);
628 p.equ(1., direction);
635 const Number p_dot_A_dot_p = p * v;
638 this->previous_alpha = alpha;
639 alpha = r_dot_preconditioner_dot_r / p_dot_A_dot_p;
642 residual_norm =
std::sqrt(std::abs(r.add_and_dot(-alpha, v, r)));
646 finalize_after_convergence(
const unsigned int)
658 template <
typename MatrixType,
typename VectorType>
659 using vmult_functions_t = decltype(std::declval<const MatrixType>().vmult(
660 std::declval<VectorType &>(),
661 std::declval<const VectorType &>(),
663 const std::function<
void(
const unsigned int,
const unsigned int)> &>(),
664 std::declval<
const std::function<
void(
const unsigned int,
665 const unsigned int)> &>()));
667 template <
typename MatrixType,
typename VectorType>
668 constexpr
bool has_vmult_functions =
669 is_supported_operation<vmult_functions_t, MatrixType, VectorType>;
674 template <
typename PreconditionerType>
675 using apply_to_subrange_t =
676 decltype(std::declval<const PreconditionerType>()
677 .apply_to_subrange(0
U, 0
U,
nullptr,
nullptr));
679 template <
typename PreconditionerType>
680 constexpr
bool has_apply_to_subrange =
681 is_supported_operation<apply_to_subrange_t, PreconditionerType>;
686 template <
typename PreconditionerType>
688 decltype(std::declval<const PreconditionerType>().
apply(0
U, 0.0));
690 template <
typename PreconditionerType>
691 constexpr
bool has_apply =
692 is_supported_operation<apply_t, PreconditionerType>;
698 template <
typename VectorType,
700 typename PreconditionerType>
701 struct IterationWorker<
705 std::enable_if_t<has_vmult_functions<MatrixType, VectorType> &&
706 (has_apply_to_subrange<PreconditionerType> ||
707 has_apply<PreconditionerType>)&&std::
708 is_same_v<VectorType,
709 LinearAlgebra::distributed::Vector<
710 typename VectorType::value_type,
713 :
public IterationWorkerBase<VectorType, MatrixType, PreconditionerType>
715 using Number =
typename VectorType::value_type;
717 Number next_r_dot_preconditioner_dot_r;
718 Number previous_beta;
720 IterationWorker(
const MatrixType &
A,
721 const PreconditionerType &preconditioner,
725 : IterationWorkerBase<VectorType, MatrixType, PreconditionerType>(
731 , next_r_dot_preconditioner_dot_r(0.)
738 do_iteration(
const unsigned int iteration_index)
740 if (iteration_index > 1)
742 previous_beta = this->beta;
743 this->beta = next_r_dot_preconditioner_dot_r /
744 this->r_dot_preconditioner_dot_r;
747 std::array<VectorizedArray<Number>, 7> vectorized_sums = {};
752 [&](
const unsigned int begin,
const unsigned int end) {
753 operation_before_loop(iteration_index,
begin,
end);
755 [&](
const unsigned int begin,
const unsigned int end) {
756 operation_after_loop(
begin,
end, vectorized_sums);
759 std::array<Number, 7> scalar_sums;
760 for (
unsigned int i = 0; i < 7; ++i)
761 scalar_sums[i] = vectorized_sums[i][0];
762 for (
unsigned int l = 1; l < VectorizedArray<Number>::size(); ++
l)
763 for (
unsigned int i = 0; i < 7; ++i)
764 scalar_sums[i] += vectorized_sums[i][
l];
768 this->r.get_mpi_communicator(),
771 this->r_dot_preconditioner_dot_r = scalar_sums[6];
773 const Number p_dot_A_dot_p = scalar_sums[0];
776 this->previous_alpha = this->alpha;
777 this->alpha = this->r_dot_preconditioner_dot_r / p_dot_A_dot_p;
781 this->residual_norm =
std::sqrt(std::abs(
783 this->alpha * (-2. * scalar_sums[2] + this->alpha * scalar_sums[1])));
785 next_r_dot_preconditioner_dot_r =
std::abs(
786 this->r_dot_preconditioner_dot_r +
787 this->alpha * (-2. * scalar_sums[4] + this->alpha * scalar_sums[5]));
792 template <
typename U =
void>
793 std::enable_if_t<has_apply<PreconditionerType>,
U>
794 operation_before_loop(
const unsigned int iteration_index,
795 const unsigned int start_range,
796 const unsigned int end_range)
const
798 Number *x = this->x.begin();
799 Number *r = this->r.begin();
800 Number *p = this->p.begin();
801 Number *v = this->v.begin();
802 const Number alpha = this->alpha;
803 const Number beta = this->beta;
805 const unsigned int end_regular =
806 start_range + (end_range - start_range) / n_lanes * n_lanes;
807 if (iteration_index == 1)
812 for (
unsigned int j = start_range; j < end_regular; j += n_lanes)
817 for (
unsigned int l = 0;
l < n_lanes; ++
l)
818 pj[
l] = this->preconditioner.apply(j +
l, rj[
l]);
823 for (
unsigned int j = end_regular; j < end_range; ++j)
825 p[j] = this->preconditioner.apply(j, r[j]);
829 else if (iteration_index % 2 == 0 && beta != Number())
831 for (
unsigned int j = start_range; j < end_regular; j += n_lanes)
839 for (
unsigned int l = 0;
l < n_lanes; ++
l)
840 prec_rj[
l] = this->preconditioner.apply(j +
l, rj[
l]);
842 pj = beta * pj + prec_rj;
847 for (
unsigned int j = end_regular; j < end_range; ++j)
849 r[j] -= alpha * v[j];
850 p[j] = beta * p[j] + this->preconditioner.apply(j, r[j]);
854 else if (iteration_index % 2 == 0 && beta == Number())
859 for (
unsigned int j = start_range; j < end_regular; j += n_lanes)
871 for (
unsigned int l = 0;
l < n_lanes; ++
l)
872 prec_rj[
l] = this->preconditioner.apply(j +
l, rj[
l]);
873 prec_rj.
store(p + j);
877 for (
unsigned int j = end_regular; j < end_range; ++j)
879 r[j] -= alpha * v[j];
880 x[j] += alpha * p[j];
881 p[j] = this->preconditioner.apply(j, r[j]);
887 const Number alpha_plus_previous_alpha_over_beta =
888 alpha + this->previous_alpha / this->previous_beta;
889 const Number previous_alpha_over_beta =
890 this->previous_alpha / this->previous_beta;
891 for (
unsigned int j = start_range; j < end_regular; j += n_lanes)
896 xj += alpha_plus_previous_alpha_over_beta * pj;
900 for (
unsigned int l = 0;
l < n_lanes; ++
l)
902 prec_rj[
l] = this->preconditioner.apply(j +
l, rj[
l]);
903 prec_vj[
l] = this->preconditioner.apply(j +
l, vj[
l]);
905 xj -= previous_alpha_over_beta * prec_rj;
909 prec_rj -= alpha * prec_vj;
910 pj = beta * pj + prec_rj;
915 for (
unsigned int j = end_regular; j < end_range; ++j)
917 x[j] += alpha_plus_previous_alpha_over_beta * p[j];
918 x[j] -= previous_alpha_over_beta *
919 this->preconditioner.apply(j, r[j]);
920 r[j] -= alpha * v[j];
921 p[j] = beta * p[j] + this->preconditioner.apply(j, r[j]);
929 template <
typename U =
void>
930 std::enable_if_t<has_apply<PreconditionerType>,
U>
931 operation_after_loop(
932 const unsigned int start_range,
933 const unsigned int end_range,
936 const Number *r = this->r.begin();
937 const Number *p = this->p.begin();
938 const Number *v = this->v.begin();
939 std::array<VectorizedArray<Number>, 7> my_sums = {};
941 const unsigned int end_regular =
942 start_range + (end_range - start_range) / n_lanes * n_lanes;
943 for (
unsigned int j = start_range; j < end_regular; j += n_lanes)
950 for (
unsigned int l = 0;
l < n_lanes; ++
l)
952 prec_vj[
l] = this->preconditioner.apply(j +
l, vj[
l]);
953 prec_rj[
l] = this->preconditioner.apply(j +
l, rj[
l]);
955 my_sums[0] += pj * vj;
956 my_sums[1] += vj * vj;
957 my_sums[2] += rj * vj;
958 my_sums[3] += rj * rj;
959 my_sums[4] += rj * prec_vj;
960 my_sums[5] += vj * prec_vj;
961 my_sums[6] += rj * prec_rj;
963 for (
unsigned int j = end_regular; j < end_range; ++j)
965 const Number prec_v = this->preconditioner.apply(j, v[j]);
966 const Number prec_r = this->preconditioner.apply(j, r[j]);
967 my_sums[0][0] += p[j] * v[j];
968 my_sums[1][0] += v[j] * v[j];
969 my_sums[2][0] += r[j] * v[j];
970 my_sums[3][0] += r[j] * r[j];
971 my_sums[4][0] += r[j] * prec_v;
972 my_sums[5][0] += v[j] * prec_v;
973 my_sums[6][0] += r[j] * prec_r;
975 for (
unsigned int i = 0; i < vectorized_sums.size(); ++i)
976 vectorized_sums[i] += my_sums[i];
981 template <
typename U =
void>
982 std::enable_if_t<has_apply<PreconditionerType>,
U>
983 finalize_after_convergence(
const unsigned int iteration_index)
985 if (iteration_index % 2 == 1 || this->beta == Number())
986 this->x.add(this->alpha, this->p);
989 using Number =
typename VectorType::value_type;
990 const unsigned int end_range = this->x.locally_owned_size();
992 Number *x = this->x.begin();
993 Number *r = this->r.begin();
994 Number *p = this->p.begin();
1000 const Number alpha_plus_previous_alpha_over_beta =
1001 this->alpha + this->previous_alpha / this->beta;
1002 const Number previous_alpha_over_beta =
1003 this->previous_alpha / this->beta;
1006 for (
unsigned int j = 0; j < end_range; ++j)
1008 x[j] += alpha_plus_previous_alpha_over_beta * p[j] -
1009 previous_alpha_over_beta *
1010 this->preconditioner.apply(j, r[j]);
1018 template <
typename U =
void>
1019 std::enable_if_t<!has_apply<PreconditionerType>,
U>
1020 operation_before_loop(
const unsigned int iteration_index,
1021 const unsigned int start_range,
1022 const unsigned int end_range)
const
1024 Number *x = this->x.begin() + start_range;
1025 Number *r = this->r.begin() + start_range;
1026 Number *p = this->p.begin() + start_range;
1027 Number *v = this->v.begin() + start_range;
1028 const Number alpha = this->alpha;
1029 const Number beta = this->beta;
1030 constexpr
unsigned int grain_size = 128;
1031 std::array<Number, grain_size> prec_r;
1032 if (iteration_index == 1)
1034 for (
unsigned int j = start_range; j < end_range; j += grain_size)
1036 const unsigned int length =
std::min(grain_size, end_range - j);
1037 this->preconditioner.apply_to_subrange(j,
1042 for (
unsigned int i = 0; i < length; ++i)
1052 else if (iteration_index % 2 == 0 && beta != Number())
1054 for (
unsigned int j = start_range; j < end_range; j += grain_size)
1056 const unsigned int length =
std::min(grain_size, end_range - j);
1058 for (
unsigned int i = 0; i < length; ++i)
1059 r[i] -= this->alpha * v[i];
1060 this->preconditioner.apply_to_subrange(j,
1065 for (
unsigned int i = 0; i < length; ++i)
1067 p[i] = this->beta * p[i] + prec_r[i];
1075 else if (iteration_index % 2 == 0 && beta == Number())
1080 for (
unsigned int j = start_range; j < end_range; j += grain_size)
1082 const unsigned int length =
std::min(grain_size, end_range - j);
1084 for (
unsigned int i = 0; i < length; ++i)
1085 r[i] -= this->alpha * v[i];
1086 this->preconditioner.apply_to_subrange(j,
1091 for (
unsigned int i = 0; i < length; ++i)
1093 x[i] += this->alpha * p[i];
1104 const Number alpha_plus_previous_alpha_over_beta =
1105 this->alpha + this->previous_alpha / this->previous_beta;
1106 const Number previous_alpha_over_beta =
1107 this->previous_alpha / this->previous_beta;
1108 for (
unsigned int j = start_range; j < end_range; j += grain_size)
1110 const unsigned int length =
std::min(grain_size, end_range - j);
1111 this->preconditioner.apply_to_subrange(j,
1116 for (
unsigned int i = 0; i < length; ++i)
1118 x[i] += alpha_plus_previous_alpha_over_beta * p[i] -
1119 previous_alpha_over_beta * prec_r[i];
1120 r[i] -= this->alpha * v[i];
1122 this->preconditioner.apply_to_subrange(j,
1127 for (
unsigned int i = 0; i < length; ++i)
1129 p[i] = this->beta * p[i] + prec_r[i];
1143 template <
typename U =
void>
1144 std::enable_if_t<!has_apply<PreconditionerType>,
U>
1145 operation_after_loop(
1146 const unsigned int start_range,
1147 const unsigned int end_range,
1150 const Number *r = this->r.begin();
1151 const Number *p = this->p.begin();
1152 const Number *v = this->v.begin();
1153 std::array<VectorizedArray<Number>, 7> my_sums = {};
1154 constexpr
unsigned int grain_size = 128;
1157 const unsigned int end_regular =
1158 start_range + (end_range - start_range) / grain_size * grain_size;
1159 std::array<Number, grain_size> prec_r;
1160 std::array<Number, grain_size> prec_v;
1161 for (
unsigned int j = start_range; j < end_regular; j += grain_size)
1163 this->preconditioner.apply_to_subrange(j,
1167 this->preconditioner.apply_to_subrange(j,
1172 for (
unsigned int i = 0; i < grain_size;
1178 prec_rj.
load(prec_r.data() + i);
1179 prec_vj.
load(prec_v.data() + i);
1181 my_sums[0] += pj * vj;
1182 my_sums[1] += vj * vj;
1183 my_sums[2] += rj * vj;
1184 my_sums[3] += rj * rj;
1185 my_sums[4] += rj * prec_vj;
1186 my_sums[5] += vj * prec_vj;
1187 my_sums[6] += rj * prec_rj;
1190 const unsigned int length = end_range - end_regular;
1192 this->preconditioner.apply_to_subrange(end_regular,
1193 end_regular + length,
1196 this->preconditioner.apply_to_subrange(end_regular,
1197 end_regular + length,
1200 for (
unsigned int j = end_regular; j < end_range; ++j)
1202 my_sums[0][0] += p[j] * v[j];
1203 my_sums[1][0] += v[j] * v[j];
1204 my_sums[2][0] += r[j] * v[j];
1205 my_sums[3][0] += r[j] * r[j];
1206 my_sums[4][0] += r[j] * prec_v[j - end_regular];
1207 my_sums[5][0] += v[j] * prec_v[j - end_regular];
1208 my_sums[6][0] += r[j] * prec_r[j - end_regular];
1210 for (
unsigned int i = 0; i < vectorized_sums.size(); ++i)
1211 vectorized_sums[i] += my_sums[i];
1217 template <
typename U =
void>
1218 std::enable_if_t<!has_apply<PreconditionerType>,
U>
1219 finalize_after_convergence(
const unsigned int iteration_index)
1221 if (iteration_index % 2 == 1 || this->beta == Number())
1222 this->x.add(this->alpha, this->p);
1225 const unsigned int end_range = this->x.locally_owned_size();
1227 Number *x = this->x.begin();
1228 Number *r = this->r.begin();
1229 Number *p = this->p.begin();
1230 const Number alpha_plus_previous_alpha_over_beta =
1231 this->alpha + this->previous_alpha / this->beta;
1232 const Number previous_alpha_over_beta =
1233 this->previous_alpha / this->beta;
1235 constexpr
unsigned int grain_size = 128;
1236 std::array<Number, grain_size> prec_r;
1237 for (
unsigned int j = 0; j < end_range; j += grain_size)
1239 const unsigned int length =
std::min(grain_size, end_range - j);
1240 this->preconditioner.apply_to_subrange(j,
1245 for (
unsigned int i = 0; i < length; ++i)
1246 x[i] += alpha_plus_previous_alpha_over_beta * p[i] -
1247 previous_alpha_over_beta * prec_r[i];
1261 template <
typename VectorType>
1262 template <
typename MatrixType,
typename PreconditionerType>
1266 const VectorType &
b,
1267 const PreconditionerType &preconditioner)
1269 using number =
typename VectorType::value_type;
1276 const bool do_eigenvalues =
1277 !condition_number_signal.empty() || !all_condition_numbers_signal.empty() ||
1278 !eigenvalues_signal.empty() || !all_eigenvalues_signal.empty();
1281 std::vector<typename VectorType::value_type>
diagonal;
1282 std::vector<typename VectorType::value_type> offdiagonal;
1284 typename VectorType::value_type eigen_beta_alpha = 0;
1288 internal::SolverCG::
1289 IterationWorker<VectorType, MatrixType, PreconditionerType>
1291 A, preconditioner, determine_beta_by_flexible_formula, this->memory, x);
1295 solver_state = this->iteration_status(0, worker.residual_norm, x);
1303 worker.do_iteration(it);
1305 print_vectors(it, x, worker.r, worker.p);
1309 this->coefficients_signal(worker.previous_alpha, worker.beta);
1314 diagonal.push_back(number(1.) / worker.previous_alpha +
1316 eigen_beta_alpha = worker.beta / worker.previous_alpha;
1317 offdiagonal.push_back(std::sqrt(worker.beta) /
1318 worker.previous_alpha);
1322 all_eigenvalues_signal,
1323 all_condition_numbers_signal);
1326 solver_state = this->iteration_status(it, worker.residual_norm, x);
1329 worker.finalize_after_convergence(it);
1334 condition_number_signal);
1342 template <
typename VectorType>
1343 boost::signals2::connection
1345 const std::function<
void(
typename VectorType::value_type,
1346 typename VectorType::value_type)> &slot)
1348 return coefficients_signal.connect(slot);
1353 template <
typename VectorType>
1354 boost::signals2::connection
1356 const std::function<
void(
double)> &slot,
1357 const bool every_iteration)
1359 if (every_iteration)
1361 return all_condition_numbers_signal.connect(slot);
1365 return condition_number_signal.connect(slot);
1371 template <
typename VectorType>
1372 boost::signals2::connection
1374 const std::function<
void(
const std::vector<double> &)> &slot,
1375 const bool every_iteration)
1377 if (every_iteration)
1379 return all_eigenvalues_signal.connect(slot);
1383 return eigenvalues_signal.connect(slot);
1389 template <
typename VectorType>
1392 const AdditionalData &)
1400 template <
typename VectorType>
1402 const AdditionalData &)
SolverCG(SolverControl &cn, VectorMemory< VectorType > &mem, const AdditionalData &data=AdditionalData())
boost::signals2::signal< void(const std::vector< double > &)> eigenvalues_signal
bool determine_beta_by_flexible_formula
virtual void print_vectors(const unsigned int step, const VectorType &x, const VectorType &r, const VectorType &d) const
boost::signals2::signal< void(typename VectorType::value_type, typename VectorType::value_type)> coefficients_signal
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::signal< void(double)> condition_number_signal
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)
boost::signals2::signal< void(double)> all_condition_numbers_signal
types::global_dof_index size_type
SolverCG(SolverControl &cn, const AdditionalData &data=AdditionalData())
boost::signals2::signal< void(const std::vector< double > &)> all_eigenvalues_signal
AdditionalData additional_data
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.
types::global_dof_index size_type
SolverFlexibleCG(SolverControl &cn, const AdditionalData &data=AdditionalData())
SolverFlexibleCG(SolverControl &cn, VectorMemory< VectorType > &mem, const AdditionalData &data=AdditionalData())
std::array< Number, 1 > eigenvalues(const SymmetricTensor< 2, 1, Number > &T)
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_NAMESPACE_CLOSE
static ::ExceptionBase & ExcDivideByZero()
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotImplemented()
#define AssertIndexRange(index, range)
#define AssertThrow(cond, exc)
void apply(const Kokkos::TeamPolicy< MemorySpace::Default::kokkos_space::execution_space >::member_type &team_member, const Kokkos::View< Number *, MemorySpace::Default::kokkos_space > shape_data, const ViewTypeIn in, ViewTypeOut out)
@ diagonal
Matrix is diagonal.
types::global_dof_index size_type
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
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 * begin(VectorType &V)
VectorType::value_type * end(VectorType &V)
T sum(const T &t, const MPI_Comm mpi_communicator)
unsigned int global_dof_index
::VectorizedArray< Number, width > abs(const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > sqrt(const ::VectorizedArray< Number, width > &)