15#ifndef dealii_trilinos_sparsity_pattern_h
16#define dealii_trilinos_sparsity_pattern_h
20#ifdef DEAL_II_WITH_TRILINOS
29# include <Epetra_FECrsGraph.h>
30# include <Epetra_Map.h>
31# include <Epetra_MpiComm.h>
117 <<
"You tried to access row " << arg1
118 <<
" of a distributed sparsity pattern, "
119 <<
" but only rows " << arg2 <<
" through " << arg3
120 <<
" are stored locally and can be accessed.");
241 <<
"Attempt to access element " << arg2 <<
" of row "
242 << arg1 <<
" which doesn't have that many elements.");
323 const std::vector<size_type> &n_entries_per_row);
369 const std::vector<size_type> &n_entries_per_row);
383 template <
typename SparsityPatternType>
385 copy_from(
const SparsityPatternType &nontrilinos_sparsity_pattern);
436 const MPI_Comm communicator = MPI_COMM_WORLD,
451 const std::vector<size_type> &n_entries_per_row);
468 const IndexSet &col_parallel_partitioning,
469 const MPI_Comm communicator = MPI_COMM_WORLD,
484 const IndexSet &col_parallel_partitioning,
486 const std::vector<size_type> &n_entries_per_row);
515 const IndexSet &col_parallel_partitioning,
517 const MPI_Comm communicator = MPI_COMM_WORLD,
537 const MPI_Comm communicator = MPI_COMM_WORLD,
553 const std::vector<size_type> &n_entries_per_row);
573 const IndexSet &col_parallel_partitioning,
574 const MPI_Comm communicator = MPI_COMM_WORLD,
604 const IndexSet &col_parallel_partitioning,
606 const MPI_Comm communicator = MPI_COMM_WORLD,
615 const IndexSet &col_parallel_partitioning,
617 const std::vector<size_type> &n_entries_per_row);
628 template <
typename SparsityPatternType>
631 const IndexSet &col_parallel_partitioning,
632 const SparsityPatternType &nontrilinos_sparsity_pattern,
633 const MPI_Comm communicator = MPI_COMM_WORLD,
634 const bool exchange_data =
false);
644 template <
typename SparsityPatternType>
647 const SparsityPatternType &nontrilinos_sparsity_pattern,
648 const MPI_Comm communicator = MPI_COMM_WORLD,
649 const bool exchange_data =
false);
689 std::pair<size_type, size_type>
768 template <
typename ForwardIterator>
771 ForwardIterator
begin,
773 const bool indices_are_sorted =
false);
778 const bool indices_are_sorted =
false)
override;
792 const Epetra_FECrsGraph &
905 print(std::ostream &out,
906 const bool write_extended_trilinos_info =
false)
const;
935 <<
"An error with error number " << arg1
936 <<
" occurred while calling a Trilinos function");
944 <<
"The entry with index <" << arg1 <<
',' << arg2
945 <<
"> does not exist.");
955 <<
"You tried to access element (" << arg1 <<
'/' << arg2
957 <<
" of a distributed matrix, but only rows in range ["
958 << arg3 <<
',' << arg4
959 <<
"] are stored locally and can be accessed.");
967 <<
"You tried to access element (" << arg1 <<
'/' << arg2
968 <<
')' <<
" of a sparse matrix, but it appears to not"
969 <<
" exist in the Trilinos sparsity pattern.");
983 std::unique_ptr<Epetra_FECrsGraph>
graph;
1008 const size_type row,
1009 const size_type index)
1014 visit_present_row();
1019 inline Accessor::size_type
1020 Accessor::row()
const
1022 Assert(a_row < sparsity_pattern->n_rows(),
1023 ExcBeyondEndOfSparsityPattern());
1029 inline Accessor::size_type
1030 Accessor::column()
const
1032 Assert(a_row < sparsity_pattern->n_rows(),
1033 ExcBeyondEndOfSparsityPattern());
1034 return (*colnum_cache)[a_index];
1039 inline Accessor::size_type
1040 Accessor::index()
const
1042 Assert(a_row < sparsity_pattern->n_rows(),
1043 ExcBeyondEndOfSparsityPattern());
1050 const size_type row,
1051 const size_type index)
1052 : accessor(sp, row,
index)
1057 inline Iterator::Iterator(
const Iterator &) =
default;
1062 Iterator::operator++()
1064 Assert(accessor.a_row < accessor.sparsity_pattern->n_rows(),
1071 if (accessor.a_index >= accessor.colnum_cache->size())
1073 accessor.a_index = 0;
1076 while (accessor.a_row < accessor.sparsity_pattern->n_rows())
1078 const auto row_length =
1079 accessor.sparsity_pattern->row_length(accessor.a_row);
1080 if (row_length == 0 ||
1081 !accessor.sparsity_pattern->row_is_stored_locally(
1088 accessor.visit_present_row();
1096 Iterator::operator++(
int)
1098 const Iterator old_state = *
this;
1105 inline const Accessor &
1106 Iterator::operator*()
const
1113 inline const Accessor *
1114 Iterator::operator->()
const
1122 Iterator::operator==(
const Iterator &other)
const
1124 return (accessor.a_row == other.accessor.a_row &&
1125 accessor.a_index == other.accessor.a_index);
1131 Iterator::operator!=(
const Iterator &other)
const
1133 return !(*
this == other);
1139 Iterator::operator<(
const Iterator &other)
const
1141 return (accessor.row() < other.accessor.row() ||
1142 (accessor.row() == other.accessor.row() &&
1143 accessor.index() < other.accessor.index()));
1153 const size_type first_valid_row = this->local_range().first;
1199 SparsityPattern::in_local_range(
const size_type index)
const
1202# ifndef DEAL_II_WITH_64BIT_INDICES
1203 begin = graph->RowMap().MinMyGID();
1204 end = graph->RowMap().MaxMyGID() + 1;
1206 begin = graph->RowMap().MinMyGID64();
1207 end = graph->RowMap().MaxMyGID64() + 1;
1219 return graph->Filled();
1240 template <
typename ForwardIterator>
1243 ForwardIterator begin,
1244 ForwardIterator end,
1263 const_cast<std::decay_t<decltype(*
begin)
> *>(&*
begin));
1270 if (row_is_stored_locally(row))
1272 graph->InsertGlobalIndices(trilinos_row_index,
n_cols, col_index_ptr);
1273 else if (nonlocal_graph.get() !=
nullptr)
1278 Assert(nonlocal_graph->RowMap().LID(
1280 ExcMessage(
"Attempted to write into off-processor matrix row "
1281 "that has not be specified as being writable upon "
1283 ierr = nonlocal_graph->InsertGlobalIndices(trilinos_row_index,
1288 ierr = graph->InsertGlobalIndices(1,
1289 &trilinos_row_index,
1298 inline const Epetra_FECrsGraph &
1299 SparsityPattern::trilinos_sparsity_pattern()
const
1307 SparsityPattern::locally_owned_domain_indices()
const
1309 return IndexSet(graph->DomainMap());
1315 SparsityPattern::locally_owned_range_indices()
const
1317 return IndexSet(graph->RangeMap());
virtual void add_entries(const ArrayView< const std::pair< size_type, size_type > > &entries)
void add_entries(const size_type row, ForwardIterator begin, ForwardIterator end, const bool indices_are_sorted=false)
bool is_compressed() const
void add(const size_type i, const size_type j)
SparsityPatternIterators::Iterator const_iterator
types::global_dof_index size_type
unsigned int row_length(const size_type row) const
std::shared_ptr< const std::vector< size_type > > colnum_cache
Accessor(const SparsityPattern *sparsity_pattern, const size_type row, const size_type index)
SparsityPattern * sparsity_pattern
bool operator==(const Iterator &) const
const Accessor & operator*() const
bool operator<(const Iterator &) const
Iterator(const SparsityPattern *sparsity_pattern, const size_type row, const size_type index)
bool operator!=(const Iterator &) const
const Accessor * operator->() const
Iterator(const Iterator &i)
IndexSet locally_owned_domain_indices() const
size_type row_length(const size_type row) const
std::unique_ptr< Epetra_FECrsGraph > graph
void add(const size_type i, const size_type j)
void print(std::ostream &out, const bool write_extended_trilinos_info=false) const
unsigned int max_entries_per_row() const
size_type bandwidth() const
void print_gnuplot(std::ostream &out) const
const_iterator end() const
MPI_Comm get_mpi_communicator() const
virtual void add_row_entries(const size_type &row, const ArrayView< const size_type > &columns, const bool indices_are_sorted=false) override
const_iterator begin(const size_type r) const
std::uint64_t n_nonzero_elements() const
std::unique_ptr< Epetra_CrsGraph > nonlocal_graph
bool exists(const size_type i, const size_type j) const
const Epetra_Map & domain_partitioner() const
std::pair< size_type, size_type > local_range() const
const Epetra_FECrsGraph & trilinos_sparsity_pattern() const
void add_entries(const size_type row, ForwardIterator begin, ForwardIterator end, const bool indices_are_sorted=false)
bool is_compressed() const
unsigned int local_size() const
void copy_from(const SparsityPattern &input_sparsity_pattern)
IndexSet locally_owned_range_indices() const
std::unique_ptr< Epetra_Map > column_space_map
const Epetra_Map & range_partitioner() const
const_iterator begin() const
const_iterator end(const size_type r) const
bool in_local_range(const size_type index) const
virtual ~SparsityPattern() override=default
std::size_t memory_consumption() const
SparsityPattern & operator=(const SparsityPattern &input_sparsity_pattern)
bool row_is_stored_locally(const size_type i) const
void reinit(const size_type m, const size_type n, const size_type n_entries_per_row=0)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcAccessToNonPresentElement(size_type arg1, size_type arg2)
#define DeclException0(Exception0)
static ::ExceptionBase & ExcBeyondEndOfSparsityPattern()
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcInvalidIndexWithinRow(size_type arg1, size_type arg2)
#define Assert(cond, exc)
static ::ExceptionBase & ExcIteratorPastEnd()
static ::ExceptionBase & ExcAccessToNonlocalRow(size_type arg1, size_type arg2, size_type arg3)
#define DeclException2(Exception2, type1, type2, outsequence)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcAccessToNonLocalElement(size_type arg1, size_type arg2, size_type arg3, size_type arg4)
static ::ExceptionBase & ExcInvalidIndex(size_type arg1, size_type arg2)
#define AssertIndexRange(index, range)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
#define DeclException1(Exception1, type1, outsequence)
static ::ExceptionBase & ExcTrilinosError(int arg1)
static ::ExceptionBase & ExcTrilinosError(int arg1)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
unsigned int global_dof_index