Reference documentation for deal.II version 9.5.0
\(\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
trilinos_solver.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2008 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16#ifndef dealii_trilinos_solver_h
17#define dealii_trilinos_solver_h
18
19
20#include <deal.II/base/config.h>
21
22#ifdef DEAL_II_WITH_TRILINOS
23
27# include <deal.II/lac/vector.h>
28
29// for AztecOO solvers
30# include <Amesos.h>
31# include <AztecOO.h>
32# include <Epetra_LinearProblem.h>
33# include <Epetra_Operator.h>
34
35// for Belos solvers
36# ifdef DEAL_II_TRILINOS_WITH_BELOS
37# include <BelosBlockCGSolMgr.hpp>
38# include <BelosBlockGmresSolMgr.hpp>
39# include <BelosEpetraAdapter.hpp>
40# include <BelosIteration.hpp>
41# include <BelosMultiVec.hpp>
42# include <BelosOperator.hpp>
43# include <BelosSolverManager.hpp>
44# endif
45
46# include <memory>
47
48
50
51namespace TrilinosWrappers
52{
53 // forward declarations
54# ifndef DOXYGEN
55 class SparseMatrix;
56 class PreconditionBase;
57# endif
58
59
79 {
80 public:
89 {
109 tfqmr
111
117 {
126 explicit AdditionalData(const bool output_solver_details = false,
127 const unsigned int gmres_restart_parameter = 30);
128
134
138 const unsigned int gmres_restart_parameter;
139 };
140
145 const AdditionalData &data = AdditionalData());
146
152 SolverControl & cn,
153 const AdditionalData &data = AdditionalData());
154
158 virtual ~SolverBase() = default;
159
165 void
166 solve(const SparseMatrix & A,
167 MPI::Vector & x,
168 const MPI::Vector & b,
169 const PreconditionBase &preconditioner);
170
178 void
179 solve(const Epetra_Operator & A,
180 MPI::Vector & x,
181 const MPI::Vector & b,
182 const PreconditionBase &preconditioner);
183
193 void
194 solve(const Epetra_Operator &A,
195 MPI::Vector & x,
196 const MPI::Vector & b,
197 const Epetra_Operator &preconditioner);
198
208 void
209 solve(const Epetra_Operator & A,
210 Epetra_MultiVector & x,
211 const Epetra_MultiVector &b,
212 const PreconditionBase & preconditioner);
213
224 void
225 solve(const Epetra_Operator & A,
226 Epetra_MultiVector & x,
227 const Epetra_MultiVector &b,
228 const Epetra_Operator & preconditioner);
229
230
231
242 void
243 solve(const SparseMatrix & A,
245 const ::Vector<double> &b,
246 const PreconditionBase & preconditioner);
247
259 void
262 const ::Vector<double> &b,
263 const PreconditionBase & preconditioner);
264
271 void
274 const ::LinearAlgebra::distributed::Vector<double> &b,
275 const PreconditionBase &preconditioner);
276
284 void
287 const ::LinearAlgebra::distributed::Vector<double> &b,
288 const PreconditionBase &preconditioner);
289
290
295 control() const;
296
301 int,
302 << "An error with error number " << arg1
303 << " occurred while calling a Trilinos function");
304
305 protected:
313
314 private:
319 template <typename Preconditioner>
320 void
321 do_solve(const Preconditioner &preconditioner);
322
326 template <typename Preconditioner>
327 void
328 set_preconditioner(AztecOO &solver, const Preconditioner &preconditioner);
329
335 std::unique_ptr<Epetra_LinearProblem> linear_problem;
336
341 std::unique_ptr<AztecOO_StatusTest> status_test;
342
347 AztecOO solver;
348
353 };
354
355
356 // provide a declaration for two explicit specializations
357 template <>
358 void
359 SolverBase::set_preconditioner(AztecOO & solver,
360 const PreconditionBase &preconditioner);
361
362 template <>
363 void
364 SolverBase::set_preconditioner(AztecOO & solver,
365 const Epetra_Operator &preconditioner);
366
367
373 class SolverCG : public SolverBase
374 {
375 public:
384 };
385
386
387
393 class SolverCGS : public SolverBase
394 {
395 public:
404 };
405
406
407
412 class SolverGMRES : public SolverBase
413 {
414 public:
423 const AdditionalData &data = AdditionalData());
424 };
425
426
427
435 {
436 public:
445 const AdditionalData &data = AdditionalData());
446 };
447
448
449
456 class SolverTFQMR : public SolverBase
457 {
458 public:
467 const AdditionalData &data = AdditionalData());
468 };
469
470
471
485 {
486 public:
492 {
496 explicit AdditionalData(const bool output_solver_details = false,
497 const std::string &solver_type = "Amesos_Klu");
498
504
523 std::string solver_type;
524 };
525
530 const AdditionalData &data = AdditionalData());
531
535 virtual ~SolverDirect() = default;
536
543 void
544 initialize(const SparseMatrix &A);
545
551 void
552 solve(MPI::Vector &x, const MPI::Vector &b);
553
559 void
561 const ::LinearAlgebra::distributed::Vector<double> &b);
562
569 void
570 solve(const SparseMatrix &A, MPI::Vector &x, const MPI::Vector &b);
571
579 void
580 solve(const SparseMatrix & A,
582 const ::Vector<double> &b);
583
590 void
593 const ::LinearAlgebra::distributed::Vector<double> &b);
594
599 control() const;
600
605 int,
606 << "An error with error number " << arg1
607 << " occurred while calling a Trilinos function");
608
609 private:
614 void
615 do_solve();
616
624
630 std::unique_ptr<Epetra_LinearProblem> linear_problem;
631
636 std::unique_ptr<Amesos_BaseSolver> solver;
637
642 };
643
644
645
646# ifdef DEAL_II_TRILINOS_WITH_BELOS
653 template <typename VectorType>
654 class SolverBelos
655 {
656 public:
661 enum SolverName
662 {
666 cg,
670 gmres,
671 } solver_name;
672
678 struct AdditionalData
679 {
683 AdditionalData(const SolverName solver_name = SolverName::cg,
684 const bool right_preconditioning = false)
685 : solver_name(solver_name)
686 , right_preconditioning(right_preconditioning)
687 {}
688
692 SolverName solver_name;
693
697 bool right_preconditioning;
698 };
699
703 SolverBelos(SolverControl & solver_control,
704 const AdditionalData & additional_data,
705 const Teuchos::RCP<Teuchos::ParameterList> &belos_parameters);
706
710 template <typename OperatorType, typename PreconditionerType>
711 void
712 solve(const OperatorType & a,
713 VectorType & x,
714 const VectorType & b,
715 const PreconditionerType &p);
716
717 private:
718 SolverControl & solver_control;
719 const AdditionalData additional_data;
720 const Teuchos::RCP<Teuchos::ParameterList> &belos_parameters;
721 };
722# endif
723
724} // namespace TrilinosWrappers
725
726
727
728# ifndef DOXYGEN
729
730# ifdef DEAL_II_TRILINOS_WITH_BELOS
731namespace TrilinosWrappers
732{
733 namespace internal
734 {
741 template <class VectorType>
742 class MultiVecWrapper
743 : public Belos::MultiVec<typename VectorType::value_type>
744 {
745 public:
749 using value_type = typename VectorType::value_type;
750
754 static bool
755 this_type_is_missing_a_specialization()
756 {
757 return false;
758 }
759
763 MultiVecWrapper(VectorType &vector)
764 {
765 this->vectors.resize(1);
766 this->vectors[0].reset(
767 &vector,
768 [](auto *) { /*Nothing to do, since vector is owned outside.*/ });
769 }
770
774 MultiVecWrapper(const VectorType &vector)
775 {
776 this->vectors.resize(1);
777 this->vectors[0].reset(
778 &const_cast<VectorType &>(vector),
779 [](auto *) { /*Nothing to do, since vector is owned outside.*/ });
780 }
781
785 virtual ~MultiVecWrapper() = default;
786
790 virtual Belos::MultiVec<value_type> *
791 Clone(const int numvecs) const
792 {
793 auto new_multi_vec = new MultiVecWrapper<VectorType>;
794
795 new_multi_vec->vectors.resize(numvecs);
796
797 for (auto &vec : new_multi_vec->vectors)
798 {
799 vec = std::make_shared<VectorType>();
800
801 AssertThrow(this->vectors.size() > 0, ExcInternalError());
802 vec->reinit(*this->vectors[0]);
803 }
804
805 return new_multi_vec;
806 }
807
811 virtual Belos::MultiVec<value_type> *
812 CloneCopy() const
813 {
815 }
816
822 virtual Belos::MultiVec<value_type> *
823 CloneCopy(const std::vector<int> &index) const
824 {
825 auto new_multi_vec = new MultiVecWrapper<VectorType>;
826
827 new_multi_vec->vectors.resize(index.size());
828
829 for (unsigned int i = 0; i < index.size(); ++i)
830 {
831 AssertThrow(static_cast<unsigned int>(index[i]) <
832 this->vectors.size(),
834
835 new_multi_vec->vectors[i] = std::make_shared<VectorType>();
836
837 AssertIndexRange(index[i], this->vectors.size());
838 *new_multi_vec->vectors[i] = *this->vectors[index[i]];
839 }
840
841 return new_multi_vec;
842 }
843
849 virtual Belos::MultiVec<value_type> *
850 CloneViewNonConst(const std::vector<int> &index)
851 {
852 auto new_multi_vec = new MultiVecWrapper<VectorType>;
853
854 new_multi_vec->vectors.resize(index.size());
855
856 for (unsigned int i = 0; i < index.size(); ++i)
857 {
858 AssertThrow(static_cast<unsigned int>(index[i]) <
859 this->vectors.size(),
861
862 new_multi_vec->vectors[i].reset(
863 this->vectors[index[i]].get(),
864 [](
865 auto
866 *) { /*Nothing to do, since we are creating only a view.*/ });
867 }
868
869 return new_multi_vec;
870 }
871
877 virtual const Belos::MultiVec<value_type> *
878 CloneView(const std::vector<int> &index) const
879 {
880 auto new_multi_vec = new MultiVecWrapper<VectorType>;
881
882 new_multi_vec->vectors.resize(index.size());
883
884 for (unsigned int i = 0; i < index.size(); ++i)
885 {
886 AssertThrow(static_cast<unsigned int>(index[i]) <
887 this->vectors.size(),
889
890 new_multi_vec->vectors[i].reset(
891 this->vectors[index[i]].get(),
892 [](
893 auto
894 *) { /*Nothing to do, since we are creating only a view.*/ });
895 }
896
897 return new_multi_vec;
898 }
899
903 virtual ptrdiff_t
904 GetGlobalLength() const
905 {
906 AssertThrow(this->vectors.size() > 0, ExcInternalError());
907
908 for (unsigned int i = 1; i < this->vectors.size(); ++i)
909 AssertDimension(this->vectors[0]->size(), this->vectors[i]->size());
910
911 return this->vectors[0]->size();
912 }
913
917 virtual int
918 GetNumberVecs() const
919 {
920 return vectors.size();
921 }
922
926 virtual void
927 MvTimesMatAddMv(const value_type alpha,
928 const Belos::MultiVec<value_type> & A_,
929 const Teuchos::SerialDenseMatrix<int, value_type> &B,
930 const value_type beta)
931 {
932 const auto &A = try_to_get_underlying_vector(A_);
933
934 const unsigned int n_rows = B.numRows();
935 const unsigned int n_cols = B.numCols();
936
937 AssertThrow(n_rows == static_cast<unsigned int>(A.GetNumberVecs()),
939 AssertThrow(n_cols == static_cast<unsigned int>(this->GetNumberVecs()),
941
942 for (unsigned int i = 0; i < n_cols; ++i)
943 (*this->vectors[i]) *= beta;
944
945 for (unsigned int i = 0; i < n_cols; ++i)
946 for (unsigned int j = 0; j < n_rows; ++j)
947 this->vectors[i]->add(alpha * B(j, i), *A.vectors[j]);
948 }
949
953 virtual void
954 MvAddMv(const value_type alpha,
955 const Belos::MultiVec<value_type> &A_,
956 const value_type beta,
957 const Belos::MultiVec<value_type> &B_)
958 {
959 const auto &A = try_to_get_underlying_vector(A_);
960 const auto &B = try_to_get_underlying_vector(B_);
961
962 AssertThrow(this->vectors.size() == A.vectors.size(),
964 AssertThrow(this->vectors.size() == B.vectors.size(),
966
967 for (unsigned int i = 0; i < this->vectors.size(); ++i)
968 {
969 this->vectors[i]->equ(alpha, *A.vectors[i]);
970 this->vectors[i]->add(beta, *B.vectors[i]);
971 }
972 }
973
977 virtual void
978 MvScale(const value_type alpha)
979 {
980 for (unsigned int i = 0; i < this->vectors.size(); ++i)
981 (*this->vectors[i]) *= alpha;
982 }
983
987 virtual void
988 MvScale(const std::vector<value_type> &alpha)
989 {
991 (void)alpha;
992 }
993
998 virtual void
999 MvTransMv(const value_type alpha,
1000 const Belos::MultiVec<value_type> & A_,
1001 Teuchos::SerialDenseMatrix<int, value_type> &B) const
1002 {
1003 const auto &A = try_to_get_underlying_vector(A_);
1004
1005 const unsigned int n_rows = B.numRows();
1006 const unsigned int n_cols = B.numCols();
1007
1008 AssertThrow(n_rows == static_cast<unsigned int>(A.GetNumberVecs()),
1010 AssertThrow(n_cols == static_cast<unsigned int>(this->GetNumberVecs()),
1012
1013 for (unsigned int i = 0; i < n_rows; ++i)
1014 for (unsigned int j = 0; j < n_cols; ++j)
1015 B(i, j) = alpha * ((*A.vectors[i]) * (*this->vectors[j]));
1016 }
1017
1022 virtual void
1023 MvDot(const Belos::MultiVec<value_type> &A_,
1024 std::vector<value_type> & b) const
1025 {
1026 const auto &A = try_to_get_underlying_vector(A_);
1027
1028 AssertThrow(this->vectors.size() == A.vectors.size(),
1030 AssertThrow(this->vectors.size() == b.size(), ExcInternalError());
1031
1032 for (unsigned int i = 0; i < this->vectors.size(); ++i)
1033 b[i] = (*this->vectors[i]) * (*A.vectors[i]);
1034 }
1035
1039 virtual void
1040 MvNorm(
1041 std::vector<typename Teuchos::ScalarTraits<value_type>::magnitudeType>
1042 & normvec,
1043 Belos::NormType type = Belos::TwoNorm) const
1044 {
1045 AssertThrow(type == Belos::TwoNorm, ExcNotImplemented());
1046 AssertThrow(this->vectors.size() == normvec.size(), ExcInternalError());
1047
1048 for (unsigned int i = 0; i < this->vectors.size(); ++i)
1049 normvec[i] = this->vectors[i]->l2_norm();
1050 }
1051
1055 virtual void
1056 SetBlock(const Belos::MultiVec<value_type> &A,
1057 const std::vector<int> & index)
1058 {
1060 (void)A;
1061 (void)index;
1062 }
1063
1067 virtual void
1068 MvRandom()
1069 {
1071 }
1072
1076 virtual void
1077 MvInit(const value_type alpha)
1078 {
1080 (void)alpha;
1081 }
1082
1086 virtual void
1087 MvPrint(std::ostream &os) const
1088 {
1090 (void)os;
1091 }
1092
1096 VectorType &
1097 genericVector()
1098 {
1099 AssertThrow(GetNumberVecs() == 1, ExcNotImplemented());
1100
1101 return *vectors[0];
1102 }
1103
1107 const VectorType &
1108 genericVector() const
1109 {
1110 AssertThrow(GetNumberVecs() == 1, ExcNotImplemented());
1111
1112 return *vectors[0];
1113 }
1114
1118 static MultiVecWrapper<VectorType> &
1119 try_to_get_underlying_vector(Belos::MultiVec<value_type> &vec_in)
1120 {
1121 auto vec = dynamic_cast<MultiVecWrapper<VectorType> *>(&vec_in);
1122
1124
1125 return *vec;
1126 }
1127
1131 const static MultiVecWrapper<VectorType> &
1132 try_to_get_underlying_vector(const Belos::MultiVec<value_type> &vec_in)
1133 {
1134 auto vec = dynamic_cast<const MultiVecWrapper<VectorType> *>(&vec_in);
1135
1137
1138 return *vec;
1139 }
1140
1141
1142# ifdef HAVE_BELOS_TSQR
1143 virtual void
1144 factorExplicit(Belos::MultiVec<value_type> &,
1145 Teuchos::SerialDenseMatrix<int, value_type> &,
1146 const bool = false)
1147 {
1148 Assert(false, ExcNotImplemented());
1149 }
1150
1151 virtual int
1152 revealRank(
1153 Teuchos::SerialDenseMatrix<int, value_type> &,
1154 const typename Teuchos::ScalarTraits<value_type>::magnitudeType &)
1155 {
1156 Assert(false, ExcNotImplemented());
1157 }
1158
1159# endif
1160
1161 private:
1165 std::vector<std::shared_ptr<VectorType>> vectors;
1166
1171 MultiVecWrapper() = default;
1172 };
1173
1180 template <class OperatorType, class VectorType>
1181 class OperatorWrapper
1182 : public Belos::Operator<typename VectorType::value_type>
1183 {
1184 public:
1188 using value_type = typename VectorType::value_type;
1189
1193 static bool
1194 this_type_is_missing_a_specialization()
1195 {
1196 return false;
1197 }
1198
1202 OperatorWrapper(const OperatorType &op)
1203 : op(op)
1204 {}
1205
1209 virtual ~OperatorWrapper() = default;
1210
1214 virtual void
1215 Apply(const Belos::MultiVec<value_type> &x,
1216 Belos::MultiVec<value_type> & y,
1217 Belos::ETrans trans = Belos::NOTRANS) const
1218 {
1219 // TODO: check for Tvmult
1220 AssertThrow(trans == Belos::NOTRANS, ExcNotImplemented());
1221
1222 op.vmult(MultiVecWrapper<VectorType>::try_to_get_underlying_vector(y)
1223 .genericVector(),
1224 MultiVecWrapper<VectorType>::try_to_get_underlying_vector(x)
1225 .genericVector());
1226 }
1227
1231 virtual bool
1232 HasApplyTranspose() const
1233 {
1234 // TODO: check for Tvmult
1235 return false;
1236 }
1237
1238 private:
1242 const OperatorType &op;
1243 };
1244
1245 } // namespace internal
1246
1247
1248
1249 template <typename VectorType>
1250 SolverBelos<VectorType>::SolverBelos(
1251 SolverControl & solver_control,
1252 const AdditionalData & additional_data,
1253 const Teuchos::RCP<Teuchos::ParameterList> &belos_parameters)
1254 : solver_control(solver_control)
1255 , additional_data(additional_data)
1256 , belos_parameters(belos_parameters)
1257 {}
1258
1259
1260
1261 template <typename VectorType>
1262 template <typename OperatorType, typename PreconditionerType>
1263 void
1264 SolverBelos<VectorType>::solve(const OperatorType & A_dealii,
1265 VectorType & x_dealii,
1266 const VectorType & b_dealii,
1267 const PreconditionerType &P_dealii)
1268 {
1269 using value_type = typename VectorType::value_type;
1270
1271 using MV = Belos::MultiVec<value_type>;
1272 using OP = Belos::Operator<value_type>;
1273
1274 Teuchos::RCP<OP> A = Teuchos::rcp(
1275 new internal::OperatorWrapper<OperatorType, VectorType>(A_dealii));
1276 Teuchos::RCP<OP> P = Teuchos::rcp(
1277 new internal::OperatorWrapper<PreconditionerType, VectorType>(P_dealii));
1278 Teuchos::RCP<MV> X =
1279 Teuchos::rcp(new internal::MultiVecWrapper<VectorType>(x_dealii));
1280 Teuchos::RCP<MV> B =
1281 Teuchos::rcp(new internal::MultiVecWrapper<VectorType>(b_dealii));
1282
1283 Teuchos::RCP<Belos::LinearProblem<value_type, MV, OP>> problem =
1284 Teuchos::rcp(new Belos::LinearProblem<value_type, MV, OP>(A, X, B));
1285
1286 if (additional_data.right_preconditioning == false)
1287 problem->setLeftPrec(P);
1288 else
1289 problem->setRightPrec(P);
1290
1291 const bool problem_set = problem->setProblem();
1292 AssertThrow(problem_set, ExcInternalError());
1293
1294 // compute initial residal
1295 VectorType r;
1296 r.reinit(x_dealii, true);
1297 A_dealii.vmult(r, x_dealii);
1298 r.sadd(-1., 1., b_dealii);
1299 const auto norm_0 = r.l2_norm();
1300
1301 if (solver_control.check(0, norm_0) != SolverControl::iterate)
1302 return;
1303
1304 double relative_tolerance_to_be_achieved =
1305 solver_control.tolerance() / norm_0;
1306 const unsigned int max_steps = solver_control.max_steps();
1307
1308 if (const auto *reduction_control =
1309 dynamic_cast<ReductionControl *>(&solver_control))
1310 relative_tolerance_to_be_achieved =
1311 std::max(relative_tolerance_to_be_achieved,
1312 reduction_control->reduction());
1313
1314 Teuchos::RCP<Teuchos::ParameterList> belos_parameters_copy(
1315 Teuchos::rcp(new Teuchos::ParameterList(*belos_parameters)));
1316
1317 belos_parameters_copy->set("Convergence Tolerance",
1318 relative_tolerance_to_be_achieved);
1319 belos_parameters_copy->set("Maximum Iterations",
1320 static_cast<int>(max_steps));
1321
1322 Teuchos::RCP<Belos::SolverManager<value_type, MV, OP>> solver;
1323
1324 if (additional_data.solver_name == SolverName::cg)
1325 solver = Teuchos::rcp(
1326 new Belos::BlockCGSolMgr<value_type, MV, OP>(problem,
1327 belos_parameters_copy));
1328 else if (additional_data.solver_name == SolverName::gmres)
1329 solver = Teuchos::rcp(
1330 new Belos::BlockGmresSolMgr<value_type, MV, OP>(problem,
1331 belos_parameters_copy));
1332 else
1334
1335 const auto flag = solver->solve();
1336
1337 solver_control.check(solver->getNumIters(), solver->achievedTol() * norm_0);
1338
1339 AssertThrow(flag == Belos::ReturnType::Converged ||
1340 ((dynamic_cast<IterationNumberControl *>(&solver_control) !=
1341 nullptr) &&
1342 (solver_control.last_step() == max_steps)),
1343 SolverControl::NoConvergence(solver_control.last_step(),
1344 solver_control.last_value()));
1345 }
1346
1347} // namespace TrilinosWrappers
1348# endif
1349
1350# endif
1351
1353
1354#endif // DEAL_II_WITH_TRILINOS
1355
1356#endif
@ iterate
Continue iteration.
virtual ~SolverBase()=default
SolverControl & control() const
void solve(const SparseMatrix &A, MPI::Vector &x, const MPI::Vector &b, const PreconditionBase &preconditioner)
const AdditionalData additional_data
std::unique_ptr< AztecOO_StatusTest > status_test
void set_preconditioner(AztecOO &solver, const Preconditioner &preconditioner)
void do_solve(const Preconditioner &preconditioner)
void solve(Epetra_Operator &A, ::LinearAlgebra::distributed::Vector< double > &x, const ::LinearAlgebra::distributed::Vector< double > &b, const PreconditionBase &preconditioner)
std::unique_ptr< Epetra_LinearProblem > linear_problem
enum TrilinosWrappers::SolverBase::SolverName solver_name
void solve(const SparseMatrix &A, ::LinearAlgebra::distributed::Vector< double > &x, const ::LinearAlgebra::distributed::Vector< double > &b, const PreconditionBase &preconditioner)
void solve(const SparseMatrix &A, ::LinearAlgebra::distributed::Vector< double > &x, const ::LinearAlgebra::distributed::Vector< double > &b)
const AdditionalData additional_data
virtual ~SolverDirect()=default
void initialize(const SparseMatrix &A)
std::unique_ptr< Amesos_BaseSolver > solver
void solve(MPI::Vector &x, const MPI::Vector &b)
std::unique_ptr< Epetra_LinearProblem > linear_problem
SolverControl & control() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcTrilinosError(int arg1)
static ::ExceptionBase & ExcTrilinosError(int arg1)
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:510
#define AssertThrow(cond, exc)
std::vector< value_type > l2_norm(const typename ::Triangulation< dim, spacedim >::cell_iterator &parent, const value_type parent_value)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)