Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
trilinos_sparsity_pattern.h
1 // ---------------------------------------------------------------------
2 // @f$Id: trilinos_sparsity_pattern.h 31019 2013-09-29 23:18:26Z bangerth @f$
3 //
4 // Copyright (C) 2008 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__trilinos_sparsity_pattern_h
18 #define __deal2__trilinos_sparsity_pattern_h
19 
20 
21 #include <deal.II/base/config.h>
22 
23 #ifdef DEAL_II_WITH_TRILINOS
24 
25 # include <deal.II/base/subscriptor.h>
26 # include <deal.II/base/index_set.h>
27 # include <deal.II/lac/exceptions.h>
28 
29 # include <vector>
30 # include <cmath>
31 # include <memory>
32 
33 # include <deal.II/base/std_cxx1x/shared_ptr.h>
34 
35 # include <Epetra_FECrsGraph.h>
36 # include <Epetra_Map.h>
37 # ifdef DEAL_II_WITH_MPI
38 # include <Epetra_MpiComm.h>
39 # include "mpi.h"
40 # else
41 # include "Epetra_SerialComm.h"
42 # endif
43 
45 
46 // forward declarations
47 class SparsityPattern;
51 
52 namespace TrilinosWrappers
53 {
54  // forward declarations
55  class SparsityPattern;
56 
57  namespace SparsityPatternIterators
58  {
59  // forward declaration
60  class Iterator;
61 
75  class Accessor
76  {
77  public:
81  typedef ::types::global_dof_index size_type;
82 
87  const size_type row,
88  const size_type index);
89 
93  Accessor (const Accessor &a);
94 
99  size_type row() const;
100 
105  size_type index() const;
106 
111  size_type column() const;
112 
116  DeclException0 (ExcBeyondEndOfSparsityPattern);
117 
121  DeclException3 (ExcAccessToNonlocalRow,
122  size_type, size_type, size_type,
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.");
127 
128  private:
133 
137  size_type a_row;
138 
142  size_type a_index;
143 
169  std_cxx1x::shared_ptr<const std::vector<size_type> > colnum_cache;
170 
179  void visit_present_row ();
180 
185  friend class Iterator;
186  };
187 
193  class Iterator
194  {
195  public:
199  typedef ::types::global_dof_index size_type;
200 
207  Iterator (const SparsityPattern *sparsity_pattern,
208  const size_type row,
209  const size_type index);
210 
214  Iterator (const Iterator &i);
215 
220 
224  Iterator operator++ (int);
225 
229  const Accessor &operator* () const;
230 
234  const Accessor *operator-> () const;
235 
241  bool operator == (const Iterator &) const;
242 
246  bool operator != (const Iterator &) const;
247 
256  bool operator < (const Iterator &) const;
257 
261  DeclException2 (ExcInvalidIndexWithinRow,
262  size_type, size_type,
263  << "Attempt to access element " << arg2
264  << " of row " << arg1
265  << " which doesn't have that many elements.");
266 
267  private:
273 
275  };
276 
277  }
278 
279 
301  {
302  public:
303 
307  typedef ::types::global_dof_index size_type;
308 
314 
324  SparsityPattern ();
325 
353  SparsityPattern (const size_type m,
354  const size_type n,
355  const size_type n_entries_per_row = 0);
356 
370  SparsityPattern (const size_type m,
371  const size_type n,
372  const std::vector<size_type> &n_entries_per_row);
373 
379  SparsityPattern (const SparsityPattern &input_sparsity_pattern);
380 
386  virtual ~SparsityPattern ();
387 
407  void
408  reinit (const size_type m,
409  const size_type n,
410  const size_type n_entries_per_row = 0);
411 
424  void
425  reinit (const size_type m,
426  const size_type n,
427  const std::vector<size_type> &n_entries_per_row);
428 
434  void
435  copy_from (const SparsityPattern &input_sparsity_pattern);
436 
444  template<typename SparsityType>
445  void
446  copy_from (const SparsityType &nontrilinos_sparsity_pattern);
447 
459  SparsityPattern &operator = (const SparsityPattern &input_sparsity_pattern);
460 
471  void clear ();
472 
490  void compress ();
492 
496 
520  SparsityPattern (const Epetra_Map &parallel_partitioning,
521  const size_type n_entries_per_row = 0);
522 
542  SparsityPattern (const Epetra_Map &parallel_partitioning,
543  const std::vector<size_type> &n_entries_per_row);
544 
573  SparsityPattern (const Epetra_Map &row_parallel_partitioning,
574  const Epetra_Map &col_parallel_partitioning,
575  const size_type n_entries_per_row = 0);
576 
597  SparsityPattern (const Epetra_Map &row_parallel_partitioning,
598  const Epetra_Map &col_parallel_partitioning,
599  const std::vector<size_type> &n_entries_per_row);
600 
626  void
627  reinit (const Epetra_Map &parallel_partitioning,
628  const size_type n_entries_per_row = 0);
629 
648  void
649  reinit (const Epetra_Map &parallel_partitioning,
650  const std::vector<size_type> &n_entries_per_row);
651 
680  void
681  reinit (const Epetra_Map &row_parallel_partitioning,
682  const Epetra_Map &col_parallel_partitioning,
683  const size_type n_entries_per_row = 0);
684 
708  void
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);
712 
728  template<typename SparsityType>
729  void
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);
734 
750  template<typename SparsityType>
751  void
752  reinit (const Epetra_Map &parallel_partitioning,
753  const SparsityType &nontrilinos_sparsity_pattern,
754  const bool exchange_data = false);
756 
760 
782  SparsityPattern (const IndexSet &parallel_partitioning,
783  const MPI_Comm &communicator = MPI_COMM_WORLD,
784  const size_type n_entries_per_row = 0);
785 
805  SparsityPattern (const IndexSet &parallel_partitioning,
806  const MPI_Comm &communicator,
807  const std::vector<size_type> &n_entries_per_row);
808 
832  SparsityPattern (const IndexSet &row_parallel_partitioning,
833  const IndexSet &col_parallel_partitioning,
834  const MPI_Comm &communicator = MPI_COMM_WORLD,
835  const size_type n_entries_per_row = 0);
836 
857  SparsityPattern (const IndexSet &row_parallel_partitioning,
858  const IndexSet &col_parallel_partitioning,
859  const MPI_Comm &communicator,
860  const std::vector<size_type> &n_entries_per_row);
861 
888  void
889  reinit (const IndexSet &parallel_partitioning,
890  const MPI_Comm &communicator = MPI_COMM_WORLD,
891  const size_type n_entries_per_row = 0);
892 
911  void
912  reinit (const IndexSet &parallel_partitioning,
913  const MPI_Comm &communicator,
914  const std::vector<size_type> &n_entries_per_row);
915 
945  void
946  reinit (const IndexSet &row_parallel_partitioning,
947  const IndexSet &col_parallel_partitioning,
948  const MPI_Comm &communicator = MPI_COMM_WORLD,
949  const size_type n_entries_per_row = 0);
950 
958  void
959  reinit (const IndexSet &row_parallel_partitioning,
960  const IndexSet &col_parallel_partitioning,
961  const MPI_Comm &communicator,
962  const std::vector<size_type> &n_entries_per_row);
963 
980  template<typename SparsityType>
981  void
982  reinit (const IndexSet &row_parallel_partitioning,
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);
987 
1003  template<typename SparsityType>
1004  void
1005  reinit (const IndexSet &parallel_partitioning,
1006  const SparsityType &nontrilinos_sparsity_pattern,
1007  const MPI_Comm &communicator = MPI_COMM_WORLD,
1008  const bool exchange_data = false);
1010 
1014 
1022  bool is_compressed () const;
1023 
1029  unsigned int max_entries_per_row () const;
1030 
1035  size_type n_rows () const;
1036 
1041  size_type n_cols () const;
1042 
1056  unsigned int local_size () const;
1057 
1072  std::pair<size_type, size_type>
1073  local_range () const;
1074 
1080  bool in_local_range (const size_type index) const;
1081 
1086  size_type n_nonzero_elements () const;
1087 
1092  size_type row_length (const size_type row) const;
1093 
1105  size_type bandwidth () const;
1106 
1113  bool empty () const;
1114 
1121  bool exists (const size_type i,
1122  const size_type j) const;
1123 
1130  std::size_t memory_consumption () const;
1131 
1133 
1141  void add (const size_type i,
1142  const size_type j);
1143 
1144 
1149  template <typename ForwardIterator>
1150  void add_entries (const size_type row,
1151  ForwardIterator begin,
1152  ForwardIterator end,
1153  const bool indices_are_sorted = false);
1155 
1159 
1166  const Epetra_FECrsGraph &trilinos_sparsity_pattern () const;
1167 
1178  const Epetra_Map &domain_partitioner () const;
1179 
1189  const Epetra_Map &range_partitioner () const;
1190 
1198  const Epetra_Map &row_partitioner () const;
1199 
1210  const Epetra_Map &col_partitioner () const;
1211 
1217  const Epetra_Comm &trilinos_communicator () const;
1219 
1223 
1228  const_iterator begin () const;
1229 
1233  const_iterator end () const;
1234 
1248  const_iterator begin (const size_type r) const;
1249 
1265  const_iterator end (const size_type r) const;
1266 
1268 
1272 
1282  void write_ascii ();
1283 
1296  void print (std::ostream &out,
1297  const bool write_extended_trilinos_info = false) const;
1298 
1324  void print_gnuplot (std::ostream &out) const;
1325 
1327 
1332  DeclException1 (ExcTrilinosError,
1333  int,
1334  << "An error with error number " << arg1
1335  << " occurred while calling a Trilinos function");
1336 
1340  DeclException2 (ExcInvalidIndex,
1341  size_type, size_type,
1342  << "The entry with index <" << arg1 << ',' << arg2
1343  << "> does not exist.");
1344 
1348  DeclException0 (ExcSourceEqualsDestination);
1349 
1353  DeclException4 (ExcAccessToNonLocalElement,
1354  size_type, size_type, size_type, size_type,
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.");
1360 
1364  DeclException2 (ExcAccessToNonPresentElement,
1365  size_type, size_type,
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.");
1371  private:
1372 
1380  std_cxx1x::shared_ptr<Epetra_Map> column_space_map;
1381 
1388 
1396  std_cxx1x::shared_ptr<Epetra_FECrsGraph> graph;
1397 
1400  };
1401 
1402 
1403 
1404 // -------------------------- inline and template functions ----------------------
1405 
1406 
1407 #ifndef DOXYGEN
1408 
1409  namespace SparsityPatternIterators
1410  {
1411 
1412  inline
1414  const size_type row,
1415  const size_type index)
1416  :
1417  sparsity_pattern(const_cast<SparsityPattern *>(sp)),
1418  a_row(row),
1419  a_index(index)
1420  {
1421  visit_present_row ();
1422  }
1423 
1424 
1425  inline
1426  Accessor::Accessor (const Accessor &a)
1427  :
1428  sparsity_pattern(a.sparsity_pattern),
1429  a_row(a.a_row),
1430  a_index(a.a_index),
1431  colnum_cache (a.colnum_cache)
1432  {}
1433 
1434 
1435  inline
1436  Accessor::size_type
1437  Accessor::row() const
1438  {
1439  Assert (a_row < sparsity_pattern->n_rows(), ExcBeyondEndOfSparsityPattern());
1440  return a_row;
1441  }
1442 
1443 
1444 
1445  inline
1446  Accessor::size_type
1447  Accessor::column() const
1448  {
1449  Assert (a_row < sparsity_pattern->n_rows(), ExcBeyondEndOfSparsityPattern());
1450  return (*colnum_cache)[a_index];
1451  }
1452 
1453 
1454 
1455  inline
1456  Accessor::size_type
1457  Accessor::index() const
1458  {
1459  Assert (a_row < sparsity_pattern->n_rows(), ExcBeyondEndOfSparsityPattern());
1460  return a_index;
1461  }
1462 
1463 
1464 
1465  inline
1466  Iterator::Iterator(const SparsityPattern *sp,
1467  const size_type row,
1468  const size_type index)
1469  :
1470  accessor(sp, row, index)
1471  {}
1472 
1473 
1474  inline
1475  Iterator::Iterator(const Iterator &i)
1476  :
1477  accessor(i.accessor)
1478  {}
1479 
1480 
1481 
1482  inline
1483  Iterator &
1484  Iterator::operator++ ()
1485  {
1486  Assert (accessor.a_row < accessor.sparsity_pattern->n_rows(),
1487  ExcIteratorPastEnd());
1488 
1489  ++accessor.a_index;
1490 
1491  // If at end of line: do one
1492  // step, then cycle until we
1493  // find a row with a nonzero
1494  // number of entries.
1495  if (accessor.a_index >= accessor.colnum_cache->size())
1496  {
1497  accessor.a_index = 0;
1498  ++accessor.a_row;
1499 
1500  while ((accessor.a_row < accessor.sparsity_pattern->n_rows())
1501  &&
1502  (accessor.sparsity_pattern->row_length(accessor.a_row) == 0))
1503  ++accessor.a_row;
1504 
1505  accessor.visit_present_row();
1506  }
1507  return *this;
1508  }
1509 
1510 
1511 
1512  inline
1513  Iterator
1514  Iterator::operator++ (int)
1515  {
1516  const Iterator old_state = *this;
1517  ++(*this);
1518  return old_state;
1519  }
1520 
1521 
1522 
1523  inline
1524  const Accessor &
1525  Iterator::operator* () const
1526  {
1527  return accessor;
1528  }
1529 
1530 
1531 
1532  inline
1533  const Accessor *
1534  Iterator::operator-> () const
1535  {
1536  return &accessor;
1537  }
1538 
1539 
1540 
1541  inline
1542  bool
1543  Iterator::operator == (const Iterator &other) const
1544  {
1545  return (accessor.a_row == other.accessor.a_row &&
1546  accessor.a_index == other.accessor.a_index);
1547  }
1548 
1549 
1550 
1551  inline
1552  bool
1553  Iterator::operator != (const Iterator &other) const
1554  {
1555  return ! (*this == other);
1556  }
1557 
1558 
1559 
1560  inline
1561  bool
1562  Iterator::operator < (const Iterator &other) const
1563  {
1564  return (accessor.row() < other.accessor.row() ||
1565  (accessor.row() == other.accessor.row() &&
1566  accessor.index() < other.accessor.index()));
1567  }
1568 
1569  }
1570 
1571 
1572 
1573  inline
1575  SparsityPattern::begin() const
1576  {
1577  return const_iterator(this, 0, 0);
1578  }
1579 
1580 
1581 
1582  inline
1584  SparsityPattern::end() const
1585  {
1586  return const_iterator(this, n_rows(), 0);
1587  }
1588 
1589 
1590 
1591  inline
1593  SparsityPattern::begin(const size_type r) const
1594  {
1595  Assert (r < n_rows(), ExcIndexRangeType<size_type>(r, 0, n_rows()));
1596  if (row_length(r) > 0)
1597  return const_iterator(this, r, 0);
1598  else
1599  return end (r);
1600  }
1601 
1602 
1603 
1604  inline
1606  SparsityPattern::end(const size_type r) const
1607  {
1608  Assert (r < n_rows(), ExcIndexRangeType<size_type>(r, 0, n_rows()));
1609 
1610  // place the iterator on the first entry
1611  // past this line, or at the end of the
1612  // matrix
1613  for (size_type i=r+1; i<n_rows(); ++i)
1614  if (row_length(i) > 0)
1615  return const_iterator(this, i, 0);
1616 
1617  // if there is no such line, then take the
1618  // end iterator of the matrix
1619  return end();
1620  }
1621 
1622 
1623 
1624  inline
1625  bool
1626  SparsityPattern::in_local_range (const size_type index) const
1627  {
1629 #ifndef DEAL_II_USE_LARGE_INDEX_TYPE
1630  begin = graph->RowMap().MinMyGID();
1631  end = graph->RowMap().MaxMyGID()+1;
1632 #else
1633  begin = graph->RowMap().MinMyGID64();
1634  end = graph->RowMap().MaxMyGID64()+1;
1635 #endif
1636 
1637  return ((index >= static_cast<size_type>(begin)) &&
1638  (index < static_cast<size_type>(end)));
1639  }
1640 
1641 
1642 
1643  inline
1644  bool
1646  {
1647  return compressed;
1648  }
1649 
1650 
1651 
1652  inline
1653  bool
1654  SparsityPattern::empty () const
1655  {
1656  return ((n_rows() == 0) && (n_cols() == 0));
1657  }
1658 
1659 
1660 
1661  inline
1662  void
1664  const size_type j)
1665  {
1666  add_entries (i, &j, &j+1);
1667  }
1668 
1669 
1670 
1671  template <typename ForwardIterator>
1672  inline
1673  void
1675  ForwardIterator begin,
1676  ForwardIterator end,
1677  const bool /*indices_are_sorted*/)
1678  {
1679  if (begin == end)
1680  return;
1681 
1682  // verify that the size of the data type Trilinos expects matches that the
1683  // iterator points to. we allow for some slippage between signed and
1684  // unsigned and only compare that they are both eiter 32 or 64 bit. to
1685  // write this test properly, not that we cannot compare the size of
1686  // '*begin' because 'begin' may be an iterator and '*begin' may be an
1687  // accessor class. consequently, we need to somehow get an actual value
1688  // from it which we can by evaluating an expression such as when
1689  // multiplying the value produced by 2
1691  sizeof((*begin)*2),
1692  ExcNotImplemented());
1693 
1694  TrilinosWrappers::types::int_type *col_index_ptr =
1695  (TrilinosWrappers::types::int_type *)(&*begin);
1696  const int n_cols = static_cast<int>(end - begin);
1697  compressed = false;
1698 
1699  const int ierr = graph->InsertGlobalIndices (1,
1701  n_cols, col_index_ptr);
1702 
1703  AssertThrow (ierr >= 0, ExcTrilinosError(ierr));
1704  }
1705 
1706 
1707 
1708  inline
1709  const Epetra_FECrsGraph &
1710  SparsityPattern::trilinos_sparsity_pattern () const
1711  {
1712  return *graph;
1713  }
1714 
1715 
1716 
1717  inline
1718  const Epetra_Map &
1719  SparsityPattern::domain_partitioner () const
1720  {
1721  return static_cast<const Epetra_Map &>(graph->DomainMap());
1722  }
1723 
1724 
1725 
1726  inline
1727  const Epetra_Map &
1728  SparsityPattern::range_partitioner () const
1729  {
1730  return static_cast<const Epetra_Map &>(graph->RangeMap());
1731  }
1732 
1733 
1734 
1735  inline
1736  const Epetra_Map &
1737  SparsityPattern::row_partitioner () const
1738  {
1739  return static_cast<const Epetra_Map &>(graph->RowMap());
1740  }
1741 
1742 
1743 
1744  inline
1745  const Epetra_Map &
1746  SparsityPattern::col_partitioner () const
1747  {
1748  return static_cast<const Epetra_Map &>(graph->ColMap());
1749  }
1750 
1751 
1752 
1753  inline
1754  const Epetra_Comm &
1755  SparsityPattern::trilinos_communicator () const
1756  {
1757  return graph->RangeMap().Comm();
1758  }
1759 
1760 
1761 
1762  inline
1763  SparsityPattern::SparsityPattern (const IndexSet &parallel_partitioning,
1764  const MPI_Comm &communicator,
1765  const size_type n_entries_per_row)
1766  :
1767  compressed (false)
1768  {
1769  Epetra_Map map = parallel_partitioning.make_trilinos_map (communicator,
1770  false);
1771  reinit (map, map, n_entries_per_row);
1772  }
1773 
1774 
1775 
1776  inline
1777  SparsityPattern::SparsityPattern (const IndexSet &parallel_partitioning,
1778  const MPI_Comm &communicator,
1779  const std::vector<size_type> &n_entries_per_row)
1780  :
1781  compressed (false)
1782  {
1783  Epetra_Map map = parallel_partitioning.make_trilinos_map (communicator,
1784  false);
1785  reinit (map, map, n_entries_per_row);
1786  }
1787 
1788 
1789 
1790  inline
1791  SparsityPattern::SparsityPattern (const IndexSet &row_parallel_partitioning,
1792  const IndexSet &col_parallel_partitioning,
1793  const MPI_Comm &communicator,
1794  const size_type n_entries_per_row)
1795  :
1796  compressed (false)
1797  {
1798  Epetra_Map row_map =
1799  row_parallel_partitioning.make_trilinos_map (communicator, false);
1800  Epetra_Map col_map =
1801  col_parallel_partitioning.make_trilinos_map (communicator, false);
1802  reinit (row_map, col_map, n_entries_per_row);
1803  }
1804 
1805 
1806 
1807  inline
1809  SparsityPattern (const IndexSet &row_parallel_partitioning,
1810  const IndexSet &col_parallel_partitioning,
1811  const MPI_Comm &communicator,
1812  const std::vector<size_type> &n_entries_per_row)
1813  :
1814  compressed (false)
1815  {
1816  Epetra_Map row_map =
1817  row_parallel_partitioning.make_trilinos_map (communicator, false);
1818  Epetra_Map col_map =
1819  col_parallel_partitioning.make_trilinos_map (communicator, false);
1820  reinit (row_map, col_map, n_entries_per_row);
1821  }
1822 
1823 
1824 
1825  inline
1826  void
1827  SparsityPattern::reinit (const IndexSet &parallel_partitioning,
1828  const MPI_Comm &communicator,
1829  const size_type n_entries_per_row)
1830  {
1831  Epetra_Map map = parallel_partitioning.make_trilinos_map (communicator,
1832  false);
1833  reinit (map, map, n_entries_per_row);
1834  }
1835 
1836 
1837 
1838  inline
1839  void SparsityPattern::reinit (const IndexSet &parallel_partitioning,
1840  const MPI_Comm &communicator,
1841  const std::vector<size_type> &n_entries_per_row)
1842  {
1843  Epetra_Map map = parallel_partitioning.make_trilinos_map (communicator,
1844  false);
1845  reinit (map, map, n_entries_per_row);
1846  }
1847 
1848 
1849 
1850  inline
1851  void SparsityPattern::reinit (const IndexSet &row_parallel_partitioning,
1852  const IndexSet &col_parallel_partitioning,
1853  const MPI_Comm &communicator,
1854  const size_type n_entries_per_row)
1855  {
1856  Epetra_Map row_map =
1857  row_parallel_partitioning.make_trilinos_map (communicator, false);
1858  Epetra_Map col_map =
1859  col_parallel_partitioning.make_trilinos_map (communicator, false);
1860  reinit (row_map, col_map, n_entries_per_row);
1861  }
1862 
1863 
1864  inline
1865  void
1866  SparsityPattern::reinit (const IndexSet &row_parallel_partitioning,
1867  const IndexSet &col_parallel_partitioning,
1868  const MPI_Comm &communicator,
1869  const std::vector<size_type> &n_entries_per_row)
1870  {
1871  Epetra_Map row_map =
1872  row_parallel_partitioning.make_trilinos_map (communicator, false);
1873  Epetra_Map col_map =
1874  col_parallel_partitioning.make_trilinos_map (communicator, false);
1875  reinit (row_map, col_map, n_entries_per_row);
1876  }
1877 
1878 
1879 
1880  template<typename SparsityType>
1881  inline
1882  void
1883  SparsityPattern::reinit (const IndexSet &row_parallel_partitioning,
1884  const IndexSet &col_parallel_partitioning,
1885  const SparsityType &nontrilinos_sparsity_pattern,
1886  const MPI_Comm &communicator,
1887  const bool exchange_data)
1888  {
1889  Epetra_Map row_map =
1890  row_parallel_partitioning.make_trilinos_map (communicator, false);
1891  Epetra_Map col_map =
1892  col_parallel_partitioning.make_trilinos_map (communicator, false);
1893  reinit (row_map, col_map, nontrilinos_sparsity_pattern, exchange_data);
1894  }
1895 
1896 
1897 
1898  template<typename SparsityType>
1899  inline
1900  void
1901  SparsityPattern::reinit (const IndexSet &parallel_partitioning,
1902  const SparsityType &nontrilinos_sparsity_pattern,
1903  const MPI_Comm &communicator,
1904  const bool exchange_data)
1905  {
1906  Epetra_Map map = parallel_partitioning.make_trilinos_map (communicator,
1907  false);
1908  reinit (map, map, nontrilinos_sparsity_pattern, exchange_data);
1909  }
1910 
1911 #endif // DOXYGEN
1912 }
1913 
1914 
1915 DEAL_II_NAMESPACE_CLOSE
1916 
1917 
1918 #endif // DEAL_II_WITH_TRILINOS
1919 
1920 
1921 /*-------------------- trilinos_sparsity_pattern.h --------------------*/
1922 
1923 #endif
1924 /*-------------------- trilinos_sparsity_pattern.h --------------------*/
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
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)
iterator begin() 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)
Definition: exceptions.h:362
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)
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
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)
Definition: exceptions.h:299
void copy_from(const SparsityPattern &input_sparsity_pattern)
size_type n_cols() const
unsigned int row_length(const size_type row) const
DeclException2(ExcInvalidIndex, size_type, size_type,<< "The entry with index <"<< arg1<< ','<< arg2<< "> does not exist.")
iterator end() const
std::pair< size_type, size_type > local_range() 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
DeclException0(ExcSourceEqualsDestination)
::ExceptionBase & ExcIteratorPastEnd()
size_type row_length(const size_type row) const
Iterator(const SparsityPattern *sparsity_pattern, const size_type row, const size_type index)
size_type n_nonzero_elements() const
size_type n_rows() const
bool empty() 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)