17 #ifndef __deal2__trilinos_sparsity_pattern_h
18 #define __deal2__trilinos_sparsity_pattern_h
21 #include <deal.II/base/config.h>
23 #ifdef DEAL_II_WITH_TRILINOS
25 # include <deal.II/base/subscriptor.h>
26 # include <deal.II/base/index_set.h>
27 # include <deal.II/lac/exceptions.h>
33 # include <deal.II/base/std_cxx1x/shared_ptr.h>
35 # include <Epetra_FECrsGraph.h>
36 # include <Epetra_Map.h>
37 # ifdef DEAL_II_WITH_MPI
38 # include <Epetra_MpiComm.h>
41 # include "Epetra_SerialComm.h"
44 DEAL_II_NAMESPACE_OPEN
52 namespace TrilinosWrappers
57 namespace SparsityPatternIterators
123 <<
"You tried to access row " << arg1
124 <<
" of a distributed sparsity pattern, "
125 <<
" but only rows " << arg2 <<
" through " << arg3
126 <<
" are stored locally and can be accessed.");
263 <<
"Attempt to access element " << arg2
264 <<
" of row " << arg1
265 <<
" which doesn't have that many elements.");
372 const std::vector<size_type> &n_entries_per_row);
427 const std::vector<size_type> &n_entries_per_row);
444 template<
typename SparsityType>
446 copy_from (
const SparsityType &nontrilinos_sparsity_pattern);
543 const std::vector<size_type> &n_entries_per_row);
574 const Epetra_Map &col_parallel_partitioning,
598 const Epetra_Map &col_parallel_partitioning,
599 const std::vector<size_type> &n_entries_per_row);
627 reinit (
const Epetra_Map ¶llel_partitioning,
649 reinit (
const Epetra_Map ¶llel_partitioning,
650 const std::vector<size_type> &n_entries_per_row);
681 reinit (
const Epetra_Map &row_parallel_partitioning,
682 const Epetra_Map &col_parallel_partitioning,
709 reinit (
const Epetra_Map &row_parallel_partitioning,
710 const Epetra_Map &col_parallel_partitioning,
711 const std::vector<size_type> &n_entries_per_row);
728 template<
typename SparsityType>
730 reinit (
const Epetra_Map &row_parallel_partitioning,
731 const Epetra_Map &col_parallel_partitioning,
732 const SparsityType &nontrilinos_sparsity_pattern,
733 const bool exchange_data =
false);
750 template<
typename SparsityType>
752 reinit (
const Epetra_Map ¶llel_partitioning,
753 const SparsityType &nontrilinos_sparsity_pattern,
754 const bool exchange_data =
false);
783 const MPI_Comm &communicator = MPI_COMM_WORLD,
806 const MPI_Comm &communicator,
807 const std::vector<size_type> &n_entries_per_row);
833 const IndexSet &col_parallel_partitioning,
834 const MPI_Comm &communicator = MPI_COMM_WORLD,
858 const IndexSet &col_parallel_partitioning,
859 const MPI_Comm &communicator,
860 const std::vector<size_type> &n_entries_per_row);
890 const MPI_Comm &communicator = MPI_COMM_WORLD,
913 const MPI_Comm &communicator,
914 const std::vector<size_type> &n_entries_per_row);
947 const IndexSet &col_parallel_partitioning,
948 const MPI_Comm &communicator = MPI_COMM_WORLD,
960 const IndexSet &col_parallel_partitioning,
961 const MPI_Comm &communicator,
962 const std::vector<size_type> &n_entries_per_row);
980 template<
typename SparsityType>
983 const IndexSet &col_parallel_partitioning,
984 const SparsityType &nontrilinos_sparsity_pattern,
985 const MPI_Comm &communicator = MPI_COMM_WORLD,
986 const bool exchange_data =
false);
1003 template<
typename SparsityType>
1006 const SparsityType &nontrilinos_sparsity_pattern,
1007 const MPI_Comm &communicator = MPI_COMM_WORLD,
1008 const bool exchange_data =
false);
1072 std::pair<size_type, size_type>
1113 bool empty ()
const;
1149 template <
typename ForwardIterator>
1151 ForwardIterator
begin,
1152 ForwardIterator
end,
1153 const bool indices_are_sorted =
false);
1296 void print (std::ostream &out,
1297 const bool write_extended_trilinos_info =
false)
const;
1334 <<
"An error with error number " << arg1
1335 <<
" occurred while calling a Trilinos function");
1342 <<
"The entry with index <" << arg1 <<
',' << arg2
1343 <<
"> does not exist.");
1355 <<
"You tried to access element (" << arg1
1356 <<
"/" << arg2 <<
")"
1357 <<
" of a distributed matrix, but only rows "
1358 << arg3 <<
" through " << arg4
1359 <<
" are stored locally and can be accessed.");
1366 <<
"You tried to access element (" << arg1
1367 <<
"/" << arg2 <<
")"
1368 <<
" of a sparse matrix, but it appears to not"
1369 <<
" exist in the Trilinos sparsity pattern.");
1396 std_cxx1x::shared_ptr<Epetra_FECrsGraph>
graph;
1409 namespace SparsityPatternIterators
1421 visit_present_row ();
1426 Accessor::Accessor (
const Accessor &a)
1428 sparsity_pattern(a.sparsity_pattern),
1431 colnum_cache (a.colnum_cache)
1437 Accessor::row()
const
1439 Assert (a_row < sparsity_pattern->n_rows(), ExcBeyondEndOfSparsityPattern());
1447 Accessor::column()
const
1449 Assert (a_row < sparsity_pattern->n_rows(), ExcBeyondEndOfSparsityPattern());
1450 return (*colnum_cache)[a_index];
1457 Accessor::index()
const
1459 Assert (a_row < sparsity_pattern->n_rows(), ExcBeyondEndOfSparsityPattern());
1470 accessor(sp, row, index)
1475 Iterator::Iterator(
const Iterator &i)
1477 accessor(i.accessor)
1484 Iterator::operator++ ()
1486 Assert (accessor.a_row < accessor.sparsity_pattern->n_rows(),
1495 if (accessor.a_index >= accessor.colnum_cache->size())
1497 accessor.a_index = 0;
1500 while ((accessor.a_row < accessor.sparsity_pattern->n_rows())
1502 (accessor.sparsity_pattern->row_length(accessor.a_row) == 0))
1505 accessor.visit_present_row();
1514 Iterator::operator++ (
int)
1516 const Iterator old_state = *
this;
1525 Iterator::operator* ()
const
1534 Iterator::operator-> ()
const
1543 Iterator::operator == (
const Iterator &other)
const
1545 return (accessor.a_row == other.accessor.a_row &&
1546 accessor.a_index == other.accessor.a_index);
1553 Iterator::operator != (
const Iterator &other)
const
1555 return ! (*
this == other);
1562 Iterator::operator < (
const Iterator &other)
const
1564 return (accessor.row() < other.accessor.row() ||
1565 (accessor.row() == other.accessor.row() &&
1566 accessor.index() < other.accessor.index()));
1626 SparsityPattern::in_local_range (
const size_type index)
const
1629 #ifndef DEAL_II_USE_LARGE_INDEX_TYPE
1630 begin = graph->RowMap().MinMyGID();
1631 end = graph->RowMap().MaxMyGID()+1;
1633 begin = graph->RowMap().MinMyGID64();
1634 end = graph->RowMap().MaxMyGID64()+1;
1637 return ((index >= static_cast<size_type>(begin)) &&
1638 (index < static_cast<size_type>(end)));
1671 template <
typename ForwardIterator>
1675 ForwardIterator begin,
1676 ForwardIterator end,
1699 const int ierr = graph->InsertGlobalIndices (1,
1701 n_cols, col_index_ptr);
1709 const Epetra_FECrsGraph &
1710 SparsityPattern::trilinos_sparsity_pattern ()
const
1719 SparsityPattern::domain_partitioner ()
const
1721 return static_cast<const Epetra_Map &
>(graph->DomainMap());
1728 SparsityPattern::range_partitioner ()
const
1730 return static_cast<const Epetra_Map &
>(graph->RangeMap());
1737 SparsityPattern::row_partitioner ()
const
1739 return static_cast<const Epetra_Map &
>(graph->RowMap());
1746 SparsityPattern::col_partitioner ()
const
1748 return static_cast<const Epetra_Map &
>(graph->ColMap());
1755 SparsityPattern::trilinos_communicator ()
const
1757 return graph->RangeMap().Comm();
1764 const MPI_Comm &communicator,
1771 reinit (map, map, n_entries_per_row);
1778 const MPI_Comm &communicator,
1779 const std::vector<size_type> &n_entries_per_row)
1785 reinit (map, map, n_entries_per_row);
1792 const IndexSet &col_parallel_partitioning,
1793 const MPI_Comm &communicator,
1798 Epetra_Map row_map =
1800 Epetra_Map col_map =
1802 reinit (row_map, col_map, n_entries_per_row);
1810 const IndexSet &col_parallel_partitioning,
1811 const MPI_Comm &communicator,
1812 const std::vector<size_type> &n_entries_per_row)
1816 Epetra_Map row_map =
1818 Epetra_Map col_map =
1820 reinit (row_map, col_map, n_entries_per_row);
1828 const MPI_Comm &communicator,
1833 reinit (map, map, n_entries_per_row);
1840 const MPI_Comm &communicator,
1841 const std::vector<size_type> &n_entries_per_row)
1845 reinit (map, map, n_entries_per_row);
1852 const IndexSet &col_parallel_partitioning,
1853 const MPI_Comm &communicator,
1856 Epetra_Map row_map =
1858 Epetra_Map col_map =
1860 reinit (row_map, col_map, n_entries_per_row);
1867 const IndexSet &col_parallel_partitioning,
1868 const MPI_Comm &communicator,
1869 const std::vector<size_type> &n_entries_per_row)
1871 Epetra_Map row_map =
1873 Epetra_Map col_map =
1875 reinit (row_map, col_map, n_entries_per_row);
1880 template<
typename SparsityType>
1884 const IndexSet &col_parallel_partitioning,
1885 const SparsityType &nontrilinos_sparsity_pattern,
1886 const MPI_Comm &communicator,
1887 const bool exchange_data)
1889 Epetra_Map row_map =
1891 Epetra_Map col_map =
1893 reinit (row_map, col_map, nontrilinos_sparsity_pattern, exchange_data);
1898 template<
typename SparsityType>
1902 const SparsityType &nontrilinos_sparsity_pattern,
1903 const MPI_Comm &communicator,
1904 const bool exchange_data)
1908 reinit (map, map, nontrilinos_sparsity_pattern, exchange_data);
1915 DEAL_II_NAMESPACE_CLOSE
1918 #endif // DEAL_II_WITH_TRILINOS
void add_entries(const size_type row, ForwardIterator begin, ForwardIterator end, const bool indices_are_sorted=false)
void print(std::ostream &out, const bool write_extended_trilinos_info=false) const
DeclException1(ExcTrilinosError, int,<< "An error with error number "<< arg1<< " occurred while calling a Trilinos function")
const_iterator end() const
std_cxx1x::shared_ptr< const std::vector< size_type > > colnum_cache
::types::global_dof_index size_type
bool operator<(const Iterator &) const
DeclException4(ExcAccessToNonLocalElement, size_type, size_type, size_type, size_type,<< "You tried to access element ("<< arg1<< "/"<< arg2<< ")"<< " of a distributed matrix, but only rows "<< arg3<< " through "<< arg4<< " are stored locally and can be accessed.")
unsigned int local_size() const
bool in_local_range(const size_type index) const
SparsityPattern & operator=(const SparsityPattern &input_sparsity_pattern)
const Accessor & operator*() const
DeclException2(ExcInvalidIndexWithinRow, size_type, size_type,<< "Attempt to access element "<< arg2<< " of row "<< arg1<< " which doesn't have that many elements.")
std_cxx1x::shared_ptr< Epetra_FECrsGraph > graph
const Epetra_Map & domain_partitioner() const
#define AssertThrow(cond, exc)
const Epetra_FECrsGraph & trilinos_sparsity_pattern() const
void add(const size_type i, const size_type j)
::types::global_dof_index size_type
void add_entries(const size_type row, ForwardIterator begin, ForwardIterator end, const bool indices_are_sorted=false)
void print_gnuplot(std::ostream &out) const
const Epetra_Map & range_partitioner() const
void reinit(const size_type m, const size_type n, const size_type n_entries_per_row=0)
const Epetra_Comm & trilinos_communicator() const
SparsityPattern * sparsity_pattern
unsigned int max_entries_per_row() const
std::size_t memory_consumption() const
const Epetra_Map & col_partitioner() const
const_iterator begin() const
bool exists(const size_type i, const size_type j) const
void add(const size_type i, const size_type j)
#define Assert(cond, exc)
void copy_from(const SparsityPattern &input_sparsity_pattern)
unsigned int row_length(const size_type row) const
DeclException2(ExcInvalidIndex, size_type, size_type,<< "The entry with index <"<< arg1<< ','<< arg2<< "> does not exist.")
const Accessor * operator->() const
::types::global_dof_index size_type
bool operator==(const Iterator &) const
std::pair< size_type, size_type > local_range() const
bool operator!=(const Iterator &) const
const Epetra_Map & row_partitioner() const
SparsityPatternIterators::Iterator const_iterator
Epetra_Map make_trilinos_map(const MPI_Comm &communicator=MPI_COMM_WORLD, const bool overlapping=false) const
SparsityPatternIterators::Iterator const_iterator
bool is_compressed() const
DeclException0(ExcSourceEqualsDestination)
size_type bandwidth() const
::ExceptionBase & ExcIteratorPastEnd()
size_type row_length(const size_type row) const
virtual ~SparsityPattern()
DeclException0(ExcBeyondEndOfSparsityPattern)
Iterator(const SparsityPattern *sparsity_pattern, const size_type row, const size_type index)
size_type n_nonzero_elements() const
::ExceptionBase & ExcNotImplemented()
bool is_compressed() const
std_cxx1x::shared_ptr< Epetra_Map > column_space_map
DeclException3(ExcAccessToNonlocalRow, size_type, size_type, size_type,<< "You tried to access row "<< arg1<< " of a distributed sparsity pattern, "<< " but only rows "<< arg2<< " through "<< arg3<< " are stored locally and can be accessed.")
Accessor(const SparsityPattern *sparsity_pattern, const size_type row, const size_type index)