16 #ifndef dealii_trilinos_sparse_matrix_h 17 # define dealii_trilinos_sparse_matrix_h 20 # include <deal.II/base/config.h> 22 # ifdef DEAL_II_WITH_TRILINOS 24 # include <deal.II/base/index_set.h> 25 # include <deal.II/base/subscriptor.h> 27 # include <deal.II/lac/exceptions.h> 28 # include <deal.II/lac/full_matrix.h> 29 # include <deal.II/lac/trilinos_epetra_vector.h> 30 # include <deal.II/lac/trilinos_vector.h> 31 # include <deal.II/lac/vector_memory.h> 32 # include <deal.II/lac/vector_operation.h> 34 # include <Epetra_Comm.h> 35 # include <Epetra_CrsGraph.h> 36 # include <Epetra_Export.h> 37 # include <Epetra_FECrsMatrix.h> 38 # include <Epetra_Map.h> 39 # include <Epetra_MultiVector.h> 40 # include <Epetra_Operator.h> 44 # include <type_traits> 46 # ifdef DEAL_II_WITH_MPI 47 # include <Epetra_MpiComm.h> 50 # include <Epetra_SerialComm.h> 53 DEAL_II_NAMESPACE_OPEN
56 template <
typename MatrixType>
59 template <
typename number>
78 template <
bool Constness>
93 <<
"You tried to access row " << arg1
94 <<
" of a distributed sparsity pattern, " 95 <<
" but only rows " << arg2 <<
" through " << arg3
96 <<
" are stored locally and can be accessed.");
197 template <
bool Constess>
236 template <
bool Other>
270 operator TrilinosScalar()
const;
276 operator=(
const TrilinosScalar n)
const;
282 operator+=(
const TrilinosScalar n)
const;
288 operator-=(
const TrilinosScalar n)
const;
294 operator*=(
const TrilinosScalar n)
const;
300 operator/=(
const TrilinosScalar n)
const;
338 friend class Reference;
355 template <
bool Constness>
379 template <
bool Other>
423 operator<(const Iterator<Constness> &)
const;
437 <<
"Attempt to access element " << arg2 <<
" of row " 438 << arg1 <<
" which doesn't have that many elements.");
446 template <
bool Other>
527 <<
"You tried to access row " << arg1
528 <<
" of a non-contiguous locally owned row set." 529 <<
" The row " << arg1
530 <<
" is not stored locally and can't be accessed.");
545 static const bool zero_addition_can_be_elided =
true;
581 const unsigned int n_max_entries_per_row);
592 const std::vector<unsigned int> &n_entries_per_row);
636 template <typename SparsityPatternType>
638 reinit(const SparsityPatternType &sparsity_pattern);
686 template <typename number>
688 reinit(const ::
SparseMatrix<number> &dealii_sparse_matrix,
689 const
double drop_tolerance = 1e-13,
690 const
bool copy_values = true,
699 reinit(const Epetra_CrsMatrix &input_matrix, const
bool copy_values = true);
721 const
size_type n_max_entries_per_row = 0);
734 const
std::vector<
unsigned int> &n_entries_per_row);
755 SparseMatrix(const Epetra_Map &row_parallel_partitioning,
756 const Epetra_Map &col_parallel_partitioning,
757 const
size_type n_max_entries_per_row = 0);
776 SparseMatrix(const Epetra_Map & row_parallel_partitioning,
777 const Epetra_Map & col_parallel_partitioning,
778 const
std::vector<
unsigned int> &n_entries_per_row);
806 template <typename SparsityPatternType>
807 DEAL_II_DEPRECATED
void 808 reinit(const Epetra_Map & parallel_partitioning,
809 const SparsityPatternType &sparsity_pattern,
810 const
bool exchange_data = false);
826 template <typename SparsityPatternType>
827 DEAL_II_DEPRECATED
void 828 reinit(const Epetra_Map & row_parallel_partitioning,
829 const Epetra_Map & col_parallel_partitioning,
830 const SparsityPatternType &sparsity_pattern,
831 const
bool exchange_data = false);
851 template <typename number>
852 DEAL_II_DEPRECATED
void 853 reinit(const Epetra_Map & parallel_partitioning,
855 const
double drop_tolerance = 1e-13,
856 const
bool copy_values = true,
874 template <typename number>
875 DEAL_II_DEPRECATED
void 876 reinit(const Epetra_Map & row_parallel_partitioning,
877 const Epetra_Map & col_parallel_partitioning,
879 const
double drop_tolerance = 1e-13,
880 const
bool copy_values = true,
900 const MPI_Comm & communicator = MPI_COMM_WORLD,
901 const
unsigned int n_max_entries_per_row = 0);
911 const MPI_Comm & communicator,
912 const
std::vector<
unsigned int> &n_entries_per_row);
929 const
IndexSet &col_parallel_partitioning,
930 const MPI_Comm &communicator = MPI_COMM_WORLD,
931 const
size_type n_max_entries_per_row = 0);
948 const
IndexSet & col_parallel_partitioning,
949 const MPI_Comm & communicator,
950 const
std::vector<
unsigned int> &n_entries_per_row);
972 template <typename SparsityPatternType>
974 reinit(const
IndexSet & parallel_partitioning,
975 const SparsityPatternType &sparsity_pattern,
976 const MPI_Comm & communicator = MPI_COMM_WORLD,
977 const
bool exchange_data = false);
991 template <typename SparsityPatternType>
993 reinit(const
IndexSet & row_parallel_partitioning,
994 const
IndexSet & col_parallel_partitioning,
995 const SparsityPatternType &sparsity_pattern,
996 const MPI_Comm & communicator = MPI_COMM_WORLD,
997 const
bool exchange_data = false);
1015 template <typename number>
1017 reinit(const
IndexSet & parallel_partitioning,
1019 const MPI_Comm & communicator = MPI_COMM_WORLD,
1020 const
double drop_tolerance = 1e-13,
1021 const
bool copy_values = true,
1037 template <typename number>
1039 reinit(const
IndexSet & row_parallel_partitioning,
1040 const
IndexSet & col_parallel_partitioning,
1042 const MPI_Comm & communicator = MPI_COMM_WORLD,
1043 const
double drop_tolerance = 1e-13,
1044 const
bool copy_values = true,
1084 local_range() const;
1091 in_local_range(const size_type
index) const;
1098 n_nonzero_elements() const;
1104 row_length(const size_type
row) const;
1113 is_compressed() const;
1121 memory_consumption() const;
1127 get_mpi_communicator() const;
1145 operator=(const
double d);
1209 set(const size_type i, const size_type j, const TrilinosScalar value);
1244 set(const
std::vector<size_type> & indices,
1245 const
FullMatrix<TrilinosScalar> &full_matrix,
1246 const
bool elide_zero_values = false);
1254 set(const
std::vector<size_type> & row_indices,
1255 const
std::vector<size_type> & col_indices,
1256 const
FullMatrix<TrilinosScalar> &full_matrix,
1257 const
bool elide_zero_values = false);
1287 set(const size_type
row,
1288 const
std::vector<size_type> & col_indices,
1289 const
std::vector<TrilinosScalar> &values,
1290 const
bool elide_zero_values = false);
1320 set(const size_type row,
1321 const size_type n_cols,
1322 const size_type * col_indices,
1323 const TrilinosScalar *values,
1324 const
bool elide_zero_values = false);
1336 add(const size_type i, const size_type j, const TrilinosScalar value);
1357 add(const
std::vector<size_type> & indices,
1358 const
FullMatrix<TrilinosScalar> &full_matrix,
1359 const
bool elide_zero_values = true);
1367 add(const
std::vector<size_type> & row_indices,
1368 const
std::vector<size_type> & col_indices,
1369 const
FullMatrix<TrilinosScalar> &full_matrix,
1370 const
bool elide_zero_values = true);
1386 add(const size_type row,
1387 const
std::vector<size_type> & col_indices,
1388 const
std::vector<TrilinosScalar> &values,
1389 const
bool elide_zero_values = true);
1405 add(const size_type row,
1406 const size_type n_cols,
1407 const size_type * col_indices,
1408 const TrilinosScalar *values,
1409 const
bool elide_zero_values = true,
1410 const
bool col_indices_are_sorted = false);
1416 operator*=(const TrilinosScalar factor);
1422 operator/=(const TrilinosScalar factor);
1467 clear_row(const size_type row, const TrilinosScalar new_diag_value = 0);
1490 clear_rows(const
std::vector<size_type> &rows,
1491 const TrilinosScalar new_diag_value = 0);
1520 operator()(const size_type i, const size_type j) const;
1539 el(const size_type i, const size_type j) const;
1548 diag_element(const size_type i) const;
1577 template <typename VectorType>
1579 vmult(VectorType &dst, const VectorType &src) const;
1603 template <typename VectorType>
1605 Tvmult(VectorType &dst, const VectorType &src) const;
1629 template <typename VectorType>
1631 vmult_add(VectorType &dst, const VectorType &src) const;
1655 template <typename VectorType>
1657 Tvmult_add(VectorType &dst, const VectorType &src) const;
1681 matrix_norm_square(const MPI::Vector &v) const;
1703 matrix_scalar_product(const MPI::Vector &u, const MPI::Vector &v) const;
1722 residual(MPI::Vector & dst,
1723 const MPI::Vector &x,
1724 const MPI::Vector &b) const;
1743 const MPI::Vector & V = MPI::Vector()) const;
1765 const MPI::Vector & V = MPI::Vector()) const;
1792 linfty_norm() const;
1799 frobenius_norm() const;
1811 const Epetra_CrsMatrix &
1812 trilinos_matrix() const;
1818 const Epetra_CrsGraph &
1819 trilinos_sparsity_pattern() const;
1830 domain_partitioner() const;
1842 range_partitioner() const;
1853 row_partitioner() const;
1866 col_partitioner() const;
1879 locally_owned_domain_indices() const;
1887 locally_owned_range_indices() const;
1965 begin(const size_type r) const;
1971 begin(const size_type r);
1983 end(const size_type r) const;
1989 end(const size_type r);
2013 print(
std::ostream &out,
2014 const
bool write_extended_trilinos_info = false) const;
2027 << "An error with error number " << arg1
2028 << " occurred while calling a Trilinos function");
2036 << "The entry with
index <" << arg1 << ',' << arg2
2037 << "> does not exist.");
2043 "You are attempting an operation on two matrices that "
2044 "are the same
object, but the operation requires that the "
2045 "two objects are in fact different.");
2060 << "You tried to access element (" << arg1 << "/" << arg2
2062 << " of a distributed matrix, but only rows " << arg3
2063 << " through " << arg4
2064 << " are stored locally and can be accessed.");
2072 << "You tried to access element (" << arg1 << "/" << arg2
2074 << " of a sparse matrix, but it appears to not"
2075 << " exist in the Trilinos sparsity pattern.");
2109 std::unique_ptr<Epetra_Map> column_space_map;
2116 std::unique_ptr<Epetra_FECrsMatrix> matrix;
2123 std::unique_ptr<Epetra_CrsMatrix> nonlocal_matrix;
2128 std::unique_ptr<Epetra_Export> nonlocal_matrix_exporter;
2141 Epetra_CombineMode last_action;
2164 check_vector_map_equality(
const Epetra_CrsMatrix & mtrx,
2165 const Epetra_MultiVector &src,
2166 const Epetra_MultiVector &dst,
2167 const bool transpose)
2169 if (transpose ==
false)
2171 Assert(src.Map().SameAs(mtrx.DomainMap()) ==
true,
2173 "Column map of matrix does not fit with vector map!"));
2174 Assert(dst.Map().SameAs(mtrx.RangeMap()) ==
true,
2175 ExcMessage(
"Row map of matrix does not fit with vector map!"));
2179 Assert(src.Map().SameAs(mtrx.RangeMap()) ==
true,
2181 "Column map of matrix does not fit with vector map!"));
2182 Assert(dst.Map().SameAs(mtrx.DomainMap()) ==
true,
2183 ExcMessage(
"Row map of matrix does not fit with vector map!"));
2191 check_vector_map_equality(
const Epetra_Operator & op,
2192 const Epetra_MultiVector &src,
2193 const Epetra_MultiVector &dst,
2194 const bool transpose)
2196 if (transpose ==
false)
2198 Assert(src.Map().SameAs(op.OperatorDomainMap()) ==
true,
2200 "Column map of operator does not fit with vector map!"));
2201 Assert(dst.Map().SameAs(op.OperatorRangeMap()) ==
true,
2203 "Row map of operator does not fit with vector map!"));
2207 Assert(src.Map().SameAs(op.OperatorRangeMap()) ==
true,
2209 "Column map of operator does not fit with vector map!"));
2210 Assert(dst.Map().SameAs(op.OperatorDomainMap()) ==
true,
2212 "Row map of operator does not fit with vector map!"));
2219 namespace LinearOperatorImplementation
2317 identity_payload()
const;
2323 null_payload()
const;
2329 transpose_payload()
const;
2347 template <
typename Solver,
typename Preconditioner>
2348 typename std::enable_if<
2349 std::is_base_of<TrilinosWrappers::SolverBase, Solver>::value &&
2351 Preconditioner>::value,
2353 inverse_payload(
Solver &,
const Preconditioner &)
const;
2372 template <
typename Solver,
typename Preconditioner>
2373 typename std::enable_if<
2374 !(std::is_base_of<TrilinosWrappers::SolverBase, Solver>::value &&
2375 std::is_base_of<TrilinosWrappers::PreconditionBase,
2376 Preconditioner>::value),
2378 inverse_payload(
Solver &,
const Preconditioner &)
const;
2393 locally_owned_domain_indices()
const;
2401 locally_owned_range_indices()
const;
2407 get_mpi_communicator()
const;
2425 std::function<void(VectorType &, const VectorType &)>
vmult;
2434 std::function<void(VectorType &, const VectorType &)>
Tvmult;
2444 std::function<void(VectorType &, const VectorType &)>
inv_vmult;
2454 std::function<void(VectorType &, const VectorType &)>
inv_Tvmult;
2470 UseTranspose()
const override;
2488 SetUseTranspose(
bool UseTranspose)
override;
2537 virtual const char *
2538 Label()
const override;
2547 virtual const Epetra_Comm &
2548 Comm()
const override;
2557 virtual const Epetra_Map &
2558 OperatorDomainMap()
const override;
2568 virtual const Epetra_Map &
2569 OperatorRangeMap()
const override;
2583 # ifdef DEAL_II_WITH_MPI 2586 Epetra_SerialComm communicator;
2610 HasNormInf()
const override;
2620 NormInf()
const override;
2656 visit_present_row();
2660 inline AccessorBase::size_type
2661 AccessorBase::row()
const 2668 inline AccessorBase::size_type
2669 AccessorBase::column()
const 2672 return (*colnum_cache)[a_index];
2676 inline AccessorBase::size_type
2677 AccessorBase::index()
const 2684 inline Accessor<true>::Accessor(MatrixType * matrix,
2685 const size_type row,
2686 const size_type index)
2687 : AccessorBase(const_cast<SparseMatrix *>(matrix), row, index)
2691 template <
bool Other>
2692 inline Accessor<true>::Accessor(
const Accessor<Other> &other)
2693 : AccessorBase(other)
2697 inline TrilinosScalar
2698 Accessor<true>::value()
const 2701 return (*value_cache)[a_index];
2705 inline Accessor<false>::Reference::Reference(
const Accessor<false> &acc)
2706 : accessor(const_cast<Accessor<false> &>(acc))
2710 inline Accessor<false>::Reference::operator TrilinosScalar()
const 2712 return (*accessor.value_cache)[accessor.a_index];
2715 inline const Accessor<false>::Reference &
2716 Accessor<false>::Reference::operator=(
const TrilinosScalar n)
const 2718 (*accessor.value_cache)[accessor.a_index] = n;
2719 accessor.matrix->set(accessor.row(),
2721 static_cast<TrilinosScalar
>(*this));
2726 inline const Accessor<false>::Reference &
2727 Accessor<false>::Reference::operator+=(
const TrilinosScalar n)
const 2729 (*accessor.value_cache)[accessor.a_index] += n;
2730 accessor.matrix->set(accessor.row(),
2732 static_cast<TrilinosScalar
>(*this));
2737 inline const Accessor<false>::Reference &
2738 Accessor<false>::Reference::operator-=(
const TrilinosScalar n)
const 2740 (*accessor.value_cache)[accessor.a_index] -= n;
2741 accessor.matrix->set(accessor.row(),
2743 static_cast<TrilinosScalar
>(*this));
2748 inline const Accessor<false>::Reference &
2749 Accessor<false>::Reference::operator*=(
const TrilinosScalar n)
const 2751 (*accessor.value_cache)[accessor.a_index] *= n;
2752 accessor.matrix->set(accessor.row(),
2754 static_cast<TrilinosScalar
>(*this));
2759 inline const Accessor<false>::Reference &
2760 Accessor<false>::Reference::operator/=(
const TrilinosScalar n)
const 2762 (*accessor.value_cache)[accessor.a_index] /= n;
2763 accessor.matrix->set(accessor.row(),
2765 static_cast<TrilinosScalar
>(*this));
2770 inline Accessor<false>::Accessor(MatrixType * matrix,
2771 const size_type row,
2772 const size_type index)
2773 : AccessorBase(matrix, row, index)
2777 inline Accessor<false>::Reference
2778 Accessor<false>::value()
const 2781 return Reference(*
this);
2786 template <
bool Constness>
2787 inline Iterator<Constness>::Iterator(MatrixType * matrix,
2788 const size_type row,
2789 const size_type index)
2790 : accessor(matrix, row, index)
2794 template <
bool Constness>
2795 template <
bool Other>
2796 inline Iterator<Constness>::Iterator(
const Iterator<Other> &other)
2797 : accessor(other.accessor)
2801 template <
bool Constness>
2802 inline Iterator<Constness> &
2803 Iterator<Constness>::operator++()
2813 if (accessor.a_index >= accessor.colnum_cache->size())
2815 accessor.a_index = 0;
2818 while ((accessor.a_row < accessor.matrix->m()) &&
2819 ((accessor.matrix->in_local_range(accessor.a_row) ==
false) ||
2820 (accessor.matrix->row_length(accessor.a_row) == 0)))
2823 accessor.visit_present_row();
2829 template <
bool Constness>
2830 inline Iterator<Constness>
2831 Iterator<Constness>::operator++(
int)
2833 const Iterator<Constness> old_state = *
this;
2840 template <
bool Constness>
2848 template <
bool Constness>
2849 inline const Accessor<Constness> *Iterator<Constness>::operator->()
const 2856 template <
bool Constness>
2858 Iterator<Constness>::operator==(
const Iterator<Constness> &other)
const 2860 return (accessor.a_row == other.accessor.a_row &&
2861 accessor.a_index == other.accessor.a_index);
2866 template <
bool Constness>
2868 Iterator<Constness>::operator!=(
const Iterator<Constness> &other)
const 2870 return !(*
this == other);
2875 template <
bool Constness>
2877 Iterator<Constness>::operator<(
const Iterator<Constness> &other)
const 2879 return (accessor.row() < other.accessor.row() ||
2880 (accessor.row() == other.accessor.row() &&
2881 accessor.index() < other.accessor.index()));
2885 template <
bool Constness>
2887 Iterator<Constness>::operator>(
const Iterator<Constness> &other)
const 2889 return (other < *
this);
2907 return const_iterator(
this, m(), 0);
2916 if (in_local_range(r) && (row_length(r) > 0))
2917 return const_iterator(
this, r, 0);
2932 for (size_type i = r + 1; i < m(); ++i)
2933 if (in_local_range(i) && (row_length(i) > 0))
2934 return const_iterator(
this, i, 0);
2954 return iterator(
this, m(), 0);
2963 if (in_local_range(r) && (row_length(r) > 0))
2964 return iterator(
this, r, 0);
2979 for (size_type i = r + 1; i < m(); ++i)
2980 if (in_local_range(i) && (row_length(i) > 0))
2981 return iterator(
this, i, 0);
2993 TrilinosWrappers::types::int_type begin, end;
2994 # ifndef DEAL_II_WITH_64BIT_INDICES 2995 begin = matrix->RowMap().MinMyGID();
2996 end = matrix->RowMap().MaxMyGID() + 1;
2998 begin = matrix->RowMap().MinMyGID64();
2999 end = matrix->RowMap().MaxMyGID64() + 1;
3002 return ((index >= static_cast<size_type>(begin)) &&
3003 (index < static_cast<size_type>(end)));
3022 const TrilinosScalar value)
3026 set(i, 1, &j, &value,
false);
3034 const bool elide_zero_values)
3036 Assert(indices.size() == values.
m(),
3040 for (size_type i = 0; i < indices.size(); ++i)
3053 const TrilinosScalar value)
3065 if (last_action == Insert)
3068 ierr = matrix->GlobalAssemble(*column_space_map,
3081 add(i, 1, &j, &value,
false);
3091 # ifndef DEAL_II_WITH_64BIT_INDICES 3092 return matrix->NumGlobalRows();
3094 return matrix->NumGlobalRows64();
3107 # ifndef DEAL_II_WITH_64BIT_INDICES 3108 return column_space_map->NumGlobalElements();
3110 return column_space_map->NumGlobalElements64();
3119 return matrix->NumMyRows();
3124 inline std::pair<SparseMatrix::size_type, SparseMatrix::size_type>
3127 size_type begin, end;
3128 # ifndef DEAL_II_WITH_64BIT_INDICES 3129 begin = matrix->RowMap().MinMyGID();
3130 end = matrix->RowMap().MaxMyGID() + 1;
3132 begin = matrix->RowMap().MinMyGID64();
3133 end = matrix->RowMap().MaxMyGID64() + 1;
3136 return std::make_pair(begin, end);
3144 # ifndef DEAL_II_WITH_64BIT_INDICES 3145 return matrix->NumGlobalNonzeros();
3147 return matrix->NumGlobalNonzeros64();
3153 template <
typename SparsityPatternType>
3156 const SparsityPatternType &sparsity_pattern,
3157 const MPI_Comm & communicator,
3158 const bool exchange_data)
3160 reinit(parallel_partitioning,
3161 parallel_partitioning,
3169 template <
typename number>
3172 const ::SparseMatrix<number> &sparse_matrix,
3173 const MPI_Comm & communicator,
3174 const double drop_tolerance,
3175 const bool copy_values,
3176 const ::SparsityPattern * use_this_sparsity)
3180 reinit(parallel_partitioning,
3181 parallel_partitioning,
3190 inline const Epetra_CrsMatrix &
3193 return static_cast<const Epetra_CrsMatrix &
>(*matrix);
3198 inline const Epetra_CrsGraph &
3201 return matrix->Graph();
3209 return IndexSet(matrix->DomainMap());
3217 return IndexSet(matrix->RangeMap());
3239 namespace LinearOperatorImplementation
3241 template <
typename Solver,
typename Preconditioner>
3242 typename std::enable_if<
3243 std::is_base_of<TrilinosWrappers::SolverBase, Solver>::value &&
3245 Preconditioner>::value,
3249 const Preconditioner &preconditioner)
const 3251 const auto &payload = *
this;
3257 return_op.inv_vmult = [payload, &solver, &preconditioner](
3262 Assert(&tril_src != &tril_dst,
3264 internal::check_vector_map_equality(payload,
3267 !payload.UseTranspose());
3268 solver.solve(payload, tril_dst, tril_src, preconditioner);
3271 return_op.inv_Tvmult = [payload, &solver, &preconditioner](
3276 Assert(&tril_src != &tril_dst,
3278 internal::check_vector_map_equality(payload,
3281 payload.UseTranspose());
3284 solver.solve(payload, tril_dst, tril_src, preconditioner);
3290 if (return_op.UseTranspose() ==
true)
3291 std::swap(return_op.inv_vmult, return_op.inv_Tvmult);
3296 template <
typename Solver,
typename Preconditioner>
3297 typename std::enable_if<
3298 !(std::is_base_of<TrilinosWrappers::SolverBase, Solver>::value &&
3299 std::is_base_of<TrilinosWrappers::PreconditionBase,
3300 Preconditioner>::value),
3309 ExcMessage(
"Payload inv_vmult disabled because of " 3310 "incompatible solver/preconditioner choice."));
3316 ExcMessage(
"Payload inv_vmult disabled because of " 3317 "incompatible solver/preconditioner choice."));
3330 DEAL_II_NAMESPACE_CLOSE
3333 # endif // DEAL_II_WITH_TRILINOS
static::ExceptionBase & ExcTrilinosError(int arg1)
#define DeclException2(Exception2, type1, type2, outsequence)
static::ExceptionBase & ExcSourceEqualsDestination()
static::ExceptionBase & ExcBeyondEndOfMatrix()
const Epetra_CrsMatrix & trilinos_matrix() const
Epetra_MultiVector VectorType
std::function< void(VectorType &, const VectorType &)> inv_Tvmult
std::function< void(VectorType &, const VectorType &)> Tvmult
static::ExceptionBase & ExcAccessToNonlocalRow(std::size_t arg1, std::size_t arg2, std::size_t arg3)
#define AssertThrow(cond, exc)
SymmetricTensor< rank_, dim, typename ProductType< Number, OtherNumber >::type > operator+(const SymmetricTensor< rank_, dim, Number > &left, const SymmetricTensor< rank_, dim, OtherNumber > &right)
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
SymmetricTensor< rank_, dim, Number > operator*(const SymmetricTensor< rank_, dim, Number > &t, const Number &factor)
const_iterator begin() const
unsigned long long int global_dof_index
bool is_compressed() const
std::shared_ptr< std::vector< size_type > > colnum_cache
::types::global_dof_index size_type
::types::global_dof_index size_type
static::ExceptionBase & ExcMessage(std::string arg1)
#define DeclException1(Exception1, type1, outsequence)
void add(const size_type i, const size_type j, const TrilinosScalar value)
#define Assert(cond, exc)
static::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define DeclExceptionMsg(Exception, defaulttext)
Epetra_MpiComm communicator
#define DeclException0(Exception0)
typename Accessor< Constness >::MatrixType MatrixType
IndexSet locally_owned_domain_indices() const
std::function< void(VectorType &, const VectorType &)> inv_vmult
std::enable_if< std::is_base_of< TrilinosWrappers::SolverBase, Solver >::value &&std::is_base_of< TrilinosWrappers::PreconditionBase, Preconditioner >::value, TrilinosPayload >::type inverse_payload(Solver &, const Preconditioner &) const
unsigned int local_size() const
SymmetricTensor< rank_, dim, Number > transpose(const SymmetricTensor< rank_, dim, Number > &t)
AccessorBase(SparseMatrix *matrix, const size_type row, const size_type index)
Epetra_Map make_trilinos_map(const MPI_Comm &communicator=MPI_COMM_WORLD, const bool overlapping=false) const
static::ExceptionBase & ExcIteratorPastEnd()
static::ExceptionBase & ExcNotQuadratic()
void swap(Vector< Number > &u, Vector< Number > &v)
Accessor< Constness > accessor
TrilinosScalar value_type
size_type n_nonzero_elements() const
void reinit(const SparsityPatternType &sparsity_pattern)
const Epetra_CrsGraph & trilinos_sparsity_pattern() const
std::function< void(VectorType &, const VectorType &)> vmult
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
const_iterator end() const
IndexSet locally_owned_range_indices() const
#define DeclException3(Exception3, type1, type2, type3, outsequence)
void set(const size_type i, const size_type j, const TrilinosScalar value)
std::pair< size_type, size_type > local_range() const
#define AssertIsFinite(number)
::types::global_dof_index size_type
bool in_local_range(const size_type index) const
static::ExceptionBase & ExcInternalError()
std::shared_ptr< std::vector< TrilinosScalar > > value_cache