16 #ifndef dealii_trilinos_sparsity_pattern_h 17 # define dealii_trilinos_sparsity_pattern_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> 29 # include <Epetra_FECrsGraph.h> 30 # include <Epetra_Map.h> 35 # ifdef DEAL_II_WITH_MPI 36 # include <Epetra_MpiComm.h> 39 # include <Epetra_SerialComm.h> 43 DEAL_II_NAMESPACE_OPEN
118 <<
"You tried to access row " << arg1
119 <<
" of a distributed sparsity pattern, " 120 <<
" but only rows " << arg2 <<
" through " << arg3
121 <<
" are stored locally and can be accessed.");
242 <<
"Attempt to access element " << arg2 <<
" of row " 243 << arg1 <<
" which doesn't have that many elements.");
325 const std::vector<size_type> &n_entries_per_row);
385 template <typename SparsityPatternType>
387 copy_from(const SparsityPatternType &nontrilinos_sparsity_pattern);
479 const Epetra_Map &col_parallel_partitioning,
497 const Epetra_Map & col_parallel_partitioning,
518 reinit(const Epetra_Map ¶llel_partitioning,
535 reinit(const Epetra_Map & parallel_partitioning,
558 reinit(const Epetra_Map &row_parallel_partitioning,
559 const Epetra_Map &col_parallel_partitioning,
577 reinit(const Epetra_Map & row_parallel_partitioning,
578 const Epetra_Map & col_parallel_partitioning,
591 template <typename SparsityPatternType>
592 DEAL_II_DEPRECATED
void 593 reinit(const Epetra_Map & row_parallel_partitioning,
594 const Epetra_Map & col_parallel_partitioning,
595 const SparsityPatternType &nontrilinos_sparsity_pattern,
596 const
bool exchange_data = false);
608 template <typename SparsityPatternType>
609 DEAL_II_DEPRECATED
void 610 reinit(const Epetra_Map & parallel_partitioning,
611 const SparsityPatternType &nontrilinos_sparsity_pattern,
612 const
bool exchange_data = false);
631 const MPI_Comm &communicator = MPI_COMM_WORLD,
645 const MPI_Comm & communicator,
663 const
IndexSet &col_parallel_partitioning,
664 const MPI_Comm &communicator = MPI_COMM_WORLD,
679 const
IndexSet & col_parallel_partitioning,
680 const MPI_Comm & communicator,
710 const
IndexSet &col_parallel_partitioning,
712 const MPI_Comm &communicator = MPI_COMM_WORLD,
731 reinit(const
IndexSet ¶llel_partitioning,
732 const MPI_Comm &communicator = MPI_COMM_WORLD,
746 reinit(const
IndexSet & parallel_partitioning,
747 const MPI_Comm & communicator,
767 reinit(const
IndexSet &row_parallel_partitioning,
768 const
IndexSet &col_parallel_partitioning,
769 const MPI_Comm &communicator = MPI_COMM_WORLD,
798 reinit(const
IndexSet &row_parallel_partitioning,
799 const
IndexSet &col_parallel_partitioning,
801 const MPI_Comm &communicator = MPI_COMM_WORLD,
809 reinit(const
IndexSet & row_parallel_partitioning,
810 const
IndexSet & col_parallel_partitioning,
811 const MPI_Comm & communicator,
823 template <typename SparsityPatternType>
825 reinit(const
IndexSet & row_parallel_partitioning,
826 const
IndexSet & col_parallel_partitioning,
827 const SparsityPatternType &nontrilinos_sparsity_pattern,
828 const MPI_Comm & communicator = MPI_COMM_WORLD,
829 const
bool exchange_data = false);
839 template <typename SparsityPatternType>
841 reinit(const
IndexSet & parallel_partitioning,
842 const SparsityPatternType &nontrilinos_sparsity_pattern,
843 const MPI_Comm & communicator = MPI_COMM_WORLD,
844 const
bool exchange_data = false);
856 is_compressed() const;
862 max_entries_per_row() const;
904 in_local_range(const size_type
index) const;
910 n_nonzero_elements() const;
921 row_length(const size_type
row) const;
944 exists(const size_type i, const size_type j) const;
951 row_is_stored_locally(const size_type i) const;
958 memory_consumption() const;
969 add(const size_type i, const size_type j);
975 template <typename ForwardIterator>
977 add_entries(const size_type
row,
978 ForwardIterator begin,
980 const
bool indices_are_sorted = false);
991 const Epetra_FECrsGraph &
992 trilinos_sparsity_pattern() const;
1004 domain_partitioner() const;
1016 range_partitioner() const;
1027 row_partitioner() const;
1040 col_partitioner() const;
1049 trilinos_communicator() const;
1055 get_mpi_communicator() const;
1069 locally_owned_domain_indices() const;
1077 locally_owned_range_indices() const;
1107 begin(const size_type r) const;
1118 end(const size_type r) const;
1142 print(
std::ostream &out,
1143 const
bool write_extended_trilinos_info = false) const;
1160 print_gnuplot(
std::ostream &out) const;
1172 << "An error with error number " << arg1
1173 << " occurred while calling a Trilinos function");
1181 << "The entry with
index <" << arg1 << ',' << arg2
1182 << "> does not exist.");
1188 ExcSourceEqualsDestination,
1189 "You are attempting an operation on two sparsity patterns that "
1190 "are the same
object, but the operation requires that the "
1191 "two objects are in fact different.");
1201 << "You tried to access element (" << arg1 << "/" << arg2
1203 << " of a distributed matrix, but only rows " << arg3
1204 << " through " << arg4
1205 << " are stored locally and can be accessed.");
1213 << "You tried to access element (" << arg1 << "/" << arg2
1215 << " of a sparse matrix, but it appears to not"
1216 << " exist in the Trilinos sparsity pattern.");
1223 std::unique_ptr<Epetra_Map> column_space_map;
1230 std::unique_ptr<Epetra_FECrsGraph> graph;
1238 std::unique_ptr<Epetra_CrsGraph> nonlocal_graph;
1255 const size_type row,
1256 const size_type
index)
1269 Assert(a_row < sparsity_pattern->n_rows(),
1279 Assert(a_row < sparsity_pattern->n_rows(),
1289 Assert(a_row < sparsity_pattern->n_rows(),
1297 const size_type row,
1298 const size_type index)
1299 : accessor(sp, row, index)
1311 Assert(accessor.a_row < accessor.sparsity_pattern->n_rows(),
1318 if (accessor.a_index >= accessor.colnum_cache->size())
1320 accessor.a_index = 0;
1323 while (accessor.a_row < accessor.sparsity_pattern->n_rows())
1325 const auto row_length =
1326 accessor.sparsity_pattern->row_length(accessor.a_row);
1327 if (row_length == 0 ||
1328 !accessor.sparsity_pattern->row_is_stored_locally(
1335 accessor.visit_present_row();
1345 const Iterator old_state = *
this;
1369 return (accessor.a_row == other.accessor.a_row &&
1370 accessor.a_index == other.accessor.a_index);
1378 return !(*
this == other);
1386 return (accessor.row() < other.accessor.row() ||
1387 (accessor.row() == other.accessor.row() &&
1388 accessor.index() < other.accessor.index()));
1398 const size_type first_valid_row = this->local_range().first;
1415 Assert(r < n_rows(), ExcIndexRangeType<size_type>(r, 0, n_rows()));
1416 if (row_length(r) > 0)
1427 Assert(r < n_rows(), ExcIndexRangeType<size_type>(r, 0, n_rows()));
1432 for (size_type i = r + 1; i < n_rows(); ++i)
1433 if (row_length(i) > 0)
1446 TrilinosWrappers::types::int_type begin, end;
1447 # ifndef DEAL_II_WITH_64BIT_INDICES 1448 begin = graph->RowMap().MinMyGID();
1449 end = graph->RowMap().MaxMyGID() + 1;
1451 begin = graph->RowMap().MinMyGID64();
1452 end = graph->RowMap().MaxMyGID64() + 1;
1455 return ((index >= static_cast<size_type>(begin)) &&
1456 (index < static_cast<size_type>(end)));
1464 return graph->Filled();
1472 return ((n_rows() == 0) && (n_cols() == 0));
1480 add_entries(i, &j, &j + 1);
1485 template <
typename ForwardIterator>
1488 ForwardIterator begin,
1489 ForwardIterator end,
1503 Assert(
sizeof(TrilinosWrappers::types::int_type) ==
sizeof((*begin) * 2),
1506 TrilinosWrappers::types::int_type *col_index_ptr =
1507 (TrilinosWrappers::types::int_type *)(&*begin);
1508 const int n_cols =
static_cast<int>(end - begin);
1511 if (row_is_stored_locally(row))
1512 ierr = graph->InsertGlobalIndices(row, n_cols, col_index_ptr);
1513 else if (nonlocal_graph.get() !=
nullptr)
1518 Assert(nonlocal_graph->RowMap().LID(
1519 static_cast<TrilinosWrappers::types::int_type>(row)) != -1,
1520 ExcMessage(
"Attempted to write into off-processor matrix row " 1521 "that has not be specified as being writable upon " 1523 ierr = nonlocal_graph->InsertGlobalIndices(row, n_cols, col_index_ptr);
1526 ierr = graph->InsertGlobalIndices(
1527 1, (TrilinosWrappers::types::int_type *)&row, n_cols, col_index_ptr);
1534 inline const Epetra_FECrsGraph &
1545 return IndexSet(static_cast<const Epetra_Map &>(graph->DomainMap()));
1553 return IndexSet(static_cast<const Epetra_Map &>(graph->RangeMap()));
1560 DEAL_II_NAMESPACE_CLOSE
1563 # endif // DEAL_II_WITH_TRILINOS static::ExceptionBase & ExcTrilinosError(int arg1)
#define DeclException2(Exception2, type1, type2, outsequence)
static::ExceptionBase & ExcBeyondEndOfSparsityPattern()
const_iterator end() const
bool operator<(const Iterator &) const
bool in_local_range(const size_type index) const
std::shared_ptr< const std::vector< size_type > > colnum_cache
::types::global_dof_index size_type
const Accessor & operator*() const
#define AssertThrow(cond, exc)
const Epetra_FECrsGraph & trilinos_sparsity_pattern() const
void add(const size_type i, const size_type j)
void add_entries(const size_type row, ForwardIterator begin, ForwardIterator end, const bool indices_are_sorted=false)
SymmetricTensor< rank_, dim, Number > operator*(const SymmetricTensor< rank_, dim, Number > &t, const Number &factor)
unsigned long long int global_dof_index
SparsityPattern * sparsity_pattern
static::ExceptionBase & ExcMessage(std::string arg1)
const_iterator begin() const
::types::global_dof_index size_type
#define DeclException1(Exception1, type1, outsequence)
IndexSet locally_owned_range_indices() const
#define Assert(cond, exc)
#define DeclExceptionMsg(Exception, defaulttext)
#define DeclException0(Exception0)
const Accessor * operator->() const
bool operator==(const Iterator &) const
bool operator!=(const Iterator &) const
IndexSet locally_owned_domain_indices() const
bool is_compressed() const
static::ExceptionBase & ExcIteratorPastEnd()
Iterator(const SparsityPattern *sparsity_pattern, const size_type row, const size_type index)
static::ExceptionBase & ExcAccessToNonlocalRow(size_type arg1, size_type arg2, size_type arg3)
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
static::ExceptionBase & ExcNotImplemented()
#define DeclException3(Exception3, type1, type2, type3, outsequence)
::types::global_dof_index size_type
Accessor(const SparsityPattern *sparsity_pattern, const size_type row, const size_type index)