Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
petsc_matrix_base.h
1 // ---------------------------------------------------------------------
2 // @f$Id: petsc_matrix_base.h 31932 2013-12-08 02:15:54Z heister @f$
3 //
4 // Copyright (C) 2004 - 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__petsc_matrix_base_h
18 #define __deal2__petsc_matrix_base_h
19 
20 
21 #include <deal.II/base/config.h>
22 
23 #ifdef DEAL_II_WITH_PETSC
24 
25 # include <deal.II/base/subscriptor.h>
26 # include <deal.II/lac/full_matrix.h>
27 # include <deal.II/lac/exceptions.h>
28 # include <deal.II/lac/vector.h>
29 
30 # include <petscmat.h>
31 # include <deal.II/base/std_cxx1x/shared_ptr.h>
32 
33 # include <vector>
34 # include <cmath>
35 
36 DEAL_II_NAMESPACE_OPEN
37 
38 template <typename Matrix> class BlockMatrixBase;
39 
40 
41 namespace PETScWrappers
42 {
43  // forward declarations
44  class VectorBase;
45  class MatrixBase;
46 
47  namespace MatrixIterators
48  {
64  {
65  private:
69  class Accessor
70  {
71  public:
76 
83  Accessor (const MatrixBase *matrix,
84  const size_type row,
85  const size_type index);
86 
90  Accessor (const Accessor &a);
91 
97  size_type row() const;
98 
104  size_type index() const;
105 
111  size_type column() const;
112 
116  PetscScalar value() const;
117 
121  DeclException0 (ExcBeyondEndOfMatrix);
125  DeclException3 (ExcAccessToNonlocalRow,
126  int, int, int,
127  << "You tried to access row " << arg1
128  << " of a distributed matrix, but only rows "
129  << arg2 << " through " << arg3
130  << " are stored locally and can be accessed.");
131 
132  private:
136  mutable MatrixBase *matrix;
137 
142 
147 
172  std_cxx1x::shared_ptr<const std::vector<size_type> > colnum_cache;
173 
178  std_cxx1x::shared_ptr<const std::vector<PetscScalar> > value_cache;
179 
187  void visit_present_row ();
188 
193  friend class const_iterator;
194  };
195 
196  public:
201 
207  const_iterator (const MatrixBase *matrix,
208  const size_type row,
209  const size_type index);
210 
215 
220 
224  const Accessor &operator* () const;
225 
229  const Accessor *operator-> () const;
230 
237  bool operator == (const const_iterator &) const;
241  bool operator != (const const_iterator &) const;
242 
252  bool operator < (const const_iterator &) const;
253 
257  DeclException2 (ExcInvalidIndexWithinRow,
258  int, int,
259  << "Attempt to access element " << arg2
260  << " of row " << arg1
261  << " which doesn't have that many elements.");
262 
263  private:
269  };
270 
271  }
272 
273 
306  class MatrixBase : public Subscriptor
307  {
308  public:
314 
319 
324  typedef PetscScalar value_type;
325 
329  MatrixBase ();
330 
335  virtual ~MatrixBase ();
336 
352  MatrixBase &
353  operator = (const value_type d);
360  void clear ();
361 
377  void set (const size_type i,
378  const size_type j,
379  const PetscScalar value);
380 
416  void set (const std::vector<size_type> &indices,
417  const FullMatrix<PetscScalar> &full_matrix,
418  const bool elide_zero_values = false);
419 
427  void set (const std::vector<size_type> &row_indices,
428  const std::vector<size_type> &col_indices,
429  const FullMatrix<PetscScalar> &full_matrix,
430  const bool elide_zero_values = false);
431 
458  void set (const size_type row,
459  const std::vector<size_type > &col_indices,
460  const std::vector<PetscScalar> &values,
461  const bool elide_zero_values = false);
462 
489  void set (const size_type row,
490  const size_type n_cols,
491  const size_type *col_indices,
492  const PetscScalar *values,
493  const bool elide_zero_values = false);
494 
510  void add (const size_type i,
511  const size_type j,
512  const PetscScalar value);
513 
551  void add (const std::vector<size_type> &indices,
552  const FullMatrix<PetscScalar> &full_matrix,
553  const bool elide_zero_values = true);
554 
562  void add (const std::vector<size_type> &row_indices,
563  const std::vector<size_type> &col_indices,
564  const FullMatrix<PetscScalar> &full_matrix,
565  const bool elide_zero_values = true);
566 
594  void add (const size_type row,
595  const std::vector<size_type> &col_indices,
596  const std::vector<PetscScalar> &values,
597  const bool elide_zero_values = true);
598 
626  void add (const size_type row,
627  const size_type n_cols,
628  const size_type *col_indices,
629  const PetscScalar *values,
630  const bool elide_zero_values = true,
631  const bool col_indices_are_sorted = false);
632 
663  void clear_row (const size_type row,
664  const PetscScalar new_diag_value = 0);
665 
679  void clear_rows (const std::vector<size_type> &rows,
680  const PetscScalar new_diag_value = 0);
681 
699  void compress (::VectorOperation::values operation);
700 
705 
722  PetscScalar operator () (const size_type i,
723  const size_type j) const;
724 
736  PetscScalar el (const size_type i,
737  const size_type j) const;
738 
754  PetscScalar diag_element (const size_type i) const;
755 
760  size_type m () const;
761 
766  size_type n () const;
767 
781  size_type local_size () const;
782 
799  std::pair<size_type, size_type>
800  local_range () const;
801 
807  bool in_local_range (const size_type index) const;
808 
815  virtual const MPI_Comm &get_mpi_communicator () const = 0;
816 
826  size_type n_nonzero_elements () const;
827 
831  size_type row_length (const size_type row) const;
832 
845  PetscReal l1_norm () const;
846 
860  PetscReal linfty_norm () const;
861 
868  PetscReal frobenius_norm () const;
869 
870 
907  PetscScalar matrix_norm_square (const VectorBase &v) const;
908 
909 
932  PetscScalar matrix_scalar_product (const VectorBase &u,
933  const VectorBase &v) const;
934 
935 
936 #if DEAL_II_PETSC_VERSION_GTE(3,1,0)
937 
943  PetscScalar trace () const;
944 #endif
945 
950  MatrixBase &operator *= (const PetscScalar factor);
951 
956  MatrixBase &operator /= (const PetscScalar factor);
957 
962  MatrixBase &add (const MatrixBase &other,
963  const PetscScalar factor);
964 
983  void vmult (VectorBase &dst,
984  const VectorBase &src) const;
985 
1006  void Tvmult (VectorBase &dst,
1007  const VectorBase &src) const;
1008 
1029  void vmult_add (VectorBase &dst,
1030  const VectorBase &src) const;
1031 
1055  void Tvmult_add (VectorBase &dst,
1056  const VectorBase &src) const;
1057 
1058 
1084  PetscScalar residual (VectorBase &dst,
1085  const VectorBase &x,
1086  const VectorBase &b) const;
1087 
1092  const_iterator begin () const;
1093 
1097  const_iterator end () const;
1098 
1111  const_iterator begin (const size_type r) const;
1112 
1124  const_iterator end (const size_type r) const;
1125 
1137  operator Mat () const;
1138 
1143  void transpose ();
1144 
1152 #if DEAL_II_PETSC_VERSION_LT(3,2,0)
1153  PetscTruth
1154 #else
1155  PetscBool
1156 #endif
1157  is_symmetric (const double tolerance = 1.e-12);
1158 
1168 #if DEAL_II_PETSC_VERSION_LT(3,2,0)
1169  PetscTruth
1170 #else
1171  PetscBool
1172 #endif
1173  is_hermitian (const double tolerance = 1.e-12);
1174 
1183  void write_ascii (const PetscViewerFormat format = PETSC_VIEWER_DEFAULT);
1184 
1189  void print (std::ostream &out,
1190  const bool alternative_output = false) const;
1191 
1196  std::size_t memory_consumption() const;
1197 
1201  DeclException1 (ExcPETScError,
1202  int,
1203  << "An error with error number " << arg1
1204  << " occurred while calling a PETSc function");
1208  DeclException0 (ExcSourceEqualsDestination);
1209 
1213  DeclException2 (ExcWrongMode,
1214  int, int,
1215  << "You tried to do a "
1216  << (arg1 == 1 ?
1217  "'set'" :
1218  (arg1 == 2 ?
1219  "'add'" : "???"))
1220  << " operation but the matrix is currently in "
1221  << (arg2 == 1 ?
1222  "'set'" :
1223  (arg2 == 2 ?
1224  "'add'" : "???"))
1225  << " mode. You first have to call 'compress()'.");
1226 
1227  protected:
1233  Mat matrix;
1234 
1256  struct LastAction
1257  {
1258  enum Values { none, insert, add };
1259  };
1260 
1265  LastAction::Values last_action;
1266 
1275  void prepare_action(const LastAction::Values new_action);
1276 
1299  void prepare_add();
1307  void prepare_set();
1308 
1309 
1310 
1311  private:
1312 
1316  MatrixBase(const MatrixBase &);
1320  MatrixBase &operator=(const MatrixBase &);
1321 
1329  std::vector<PetscInt> column_indices;
1330 
1338  std::vector<PetscScalar> column_values;
1339 
1340 
1346  template <class> friend class ::BlockMatrixBase;
1347  };
1348 
1349 
1350 
1351 #ifndef DOXYGEN
1352 // -------------------------- inline and template functions ----------------------
1353 
1354 
1355  namespace MatrixIterators
1356  {
1357 
1358  inline
1360  Accessor (const MatrixBase *matrix,
1361  const size_type row,
1362  const size_type index)
1363  :
1364  matrix(const_cast<MatrixBase *>(matrix)),
1365  a_row(row),
1366  a_index(index)
1367  {
1368  visit_present_row ();
1369  }
1370 
1371 
1372  inline
1374  Accessor (const Accessor &a)
1375  :
1376  matrix(a.matrix),
1377  a_row(a.a_row),
1378  a_index(a.a_index),
1379  colnum_cache (a.colnum_cache),
1380  value_cache (a.value_cache)
1381  {}
1382 
1383 
1384  inline
1385  const_iterator::Accessor::size_type
1386  const_iterator::Accessor::row() const
1387  {
1388  Assert (a_row < matrix->m(), ExcBeyondEndOfMatrix());
1389  return a_row;
1390  }
1391 
1392 
1393  inline
1394  const_iterator::Accessor::size_type
1395  const_iterator::Accessor::column() const
1396  {
1397  Assert (a_row < matrix->m(), ExcBeyondEndOfMatrix());
1398  return (*colnum_cache)[a_index];
1399  }
1400 
1401 
1402  inline
1403  const_iterator::Accessor::size_type
1404  const_iterator::Accessor::index() const
1405  {
1406  Assert (a_row < matrix->m(), ExcBeyondEndOfMatrix());
1407  return a_index;
1408  }
1409 
1410 
1411  inline
1412  PetscScalar
1413  const_iterator::Accessor::value() const
1414  {
1415  Assert (a_row < matrix->m(), ExcBeyondEndOfMatrix());
1416  return (*value_cache)[a_index];
1417  }
1418 
1419 
1420  inline
1421  const_iterator::
1422  const_iterator(const MatrixBase *matrix,
1423  const size_type row,
1424  const size_type index)
1425  :
1426  accessor(matrix, row, index)
1427  {}
1428 
1429 
1430 
1431  inline
1432  const_iterator &
1433  const_iterator::operator++ ()
1434  {
1435  Assert (accessor.a_row < accessor.matrix->m(), ExcIteratorPastEnd());
1436 
1437  ++accessor.a_index;
1438 
1439  // if at end of line: do one step, then
1440  // cycle until we find a row with a
1441  // nonzero number of entries
1442  if (accessor.a_index >= accessor.colnum_cache->size())
1443  {
1444  accessor.a_index = 0;
1445  ++accessor.a_row;
1446 
1447  while ((accessor.a_row < accessor.matrix->m())
1448  &&
1449  (accessor.matrix->row_length(accessor.a_row) == 0))
1450  ++accessor.a_row;
1451 
1452  accessor.visit_present_row();
1453  }
1454  return *this;
1455  }
1456 
1457 
1458  inline
1459  const_iterator
1460  const_iterator::operator++ (int)
1461  {
1462  const const_iterator old_state = *this;
1463  ++(*this);
1464  return old_state;
1465  }
1466 
1467 
1468  inline
1469  const const_iterator::Accessor &
1470  const_iterator::operator* () const
1471  {
1472  return accessor;
1473  }
1474 
1475 
1476  inline
1477  const const_iterator::Accessor *
1478  const_iterator::operator-> () const
1479  {
1480  return &accessor;
1481  }
1482 
1483 
1484  inline
1485  bool
1486  const_iterator::
1487  operator == (const const_iterator &other) const
1488  {
1489  return (accessor.a_row == other.accessor.a_row &&
1490  accessor.a_index == other.accessor.a_index);
1491  }
1492 
1493 
1494  inline
1495  bool
1496  const_iterator::
1497  operator != (const const_iterator &other) const
1498  {
1499  return ! (*this == other);
1500  }
1501 
1502 
1503  inline
1504  bool
1505  const_iterator::
1506  operator < (const const_iterator &other) const
1507  {
1508  return (accessor.row() < other.accessor.row() ||
1509  (accessor.row() == other.accessor.row() &&
1510  accessor.index() < other.accessor.index()));
1511  }
1512 
1513  }
1514 
1515 
1516 
1517  // Inline the set() and add()
1518  // functions, since they will be
1519  // called frequently, and the
1520  // compiler can optimize away
1521  // some unnecessary loops when
1522  // the sizes are given at
1523  // compile time.
1524  inline
1525  void
1526  MatrixBase::set (const size_type i,
1527  const size_type j,
1528  const PetscScalar value)
1529  {
1531 
1532  set (i, 1, &j, &value, false);
1533  }
1534 
1535 
1536 
1537  inline
1538  void
1539  MatrixBase::set (const std::vector<size_type> &indices,
1540  const FullMatrix<PetscScalar> &values,
1541  const bool elide_zero_values)
1542  {
1543  Assert (indices.size() == values.m(),
1544  ExcDimensionMismatch(indices.size(), values.m()));
1545  Assert (values.m() == values.n(), ExcNotQuadratic());
1546 
1547  for (size_type i=0; i<indices.size(); ++i)
1548  set (indices[i], indices.size(), &indices[0], &values(i,0),
1549  elide_zero_values);
1550  }
1551 
1552 
1553 
1554  inline
1555  void
1556  MatrixBase::set (const std::vector<size_type> &row_indices,
1557  const std::vector<size_type> &col_indices,
1558  const FullMatrix<PetscScalar> &values,
1559  const bool elide_zero_values)
1560  {
1561  Assert (row_indices.size() == values.m(),
1562  ExcDimensionMismatch(row_indices.size(), values.m()));
1563  Assert (col_indices.size() == values.n(),
1564  ExcDimensionMismatch(col_indices.size(), values.n()));
1565 
1566  for (size_type i=0; i<row_indices.size(); ++i)
1567  set (row_indices[i], col_indices.size(), &col_indices[0], &values(i,0),
1568  elide_zero_values);
1569  }
1570 
1571 
1572 
1573  inline
1574  void
1575  MatrixBase::set (const size_type row,
1576  const std::vector<size_type> &col_indices,
1577  const std::vector<PetscScalar> &values,
1578  const bool elide_zero_values)
1579  {
1580  Assert (col_indices.size() == values.size(),
1581  ExcDimensionMismatch(col_indices.size(), values.size()));
1582 
1583  set (row, col_indices.size(), &col_indices[0], &values[0],
1584  elide_zero_values);
1585  }
1586 
1587 
1588 
1589  inline
1590  void
1591  MatrixBase::set (const size_type row,
1592  const size_type n_cols,
1593  const size_type *col_indices,
1594  const PetscScalar *values,
1595  const bool elide_zero_values)
1596  {
1597  prepare_action(LastAction::insert);
1598 
1599  const PetscInt petsc_i = row;
1600  PetscInt *col_index_ptr;
1601 
1602  PetscScalar const *col_value_ptr;
1603  int n_columns;
1604 
1605  // If we don't elide zeros, the pointers
1606  // are already available...
1607 #ifndef PETSC_USE_64BIT_INDICES
1608  if (elide_zero_values == false)
1609  {
1610  col_index_ptr = (int *)col_indices;
1611  col_value_ptr = values;
1612  n_columns = n_cols;
1613  }
1614  else
1615 #endif
1616  {
1617  // Otherwise, extract nonzero values in
1618  // each row and get the respective index.
1619  if (column_indices.size() < n_cols)
1620  {
1621  column_indices.resize(n_cols);
1622  column_values.resize(n_cols);
1623  }
1624 
1625  n_columns = 0;
1626  for (size_type j=0; j<n_cols; ++j)
1627  {
1628  const PetscScalar value = values[j];
1630  if (value != PetscScalar())
1631  {
1632  column_indices[n_columns] = col_indices[j];
1633  column_values[n_columns] = value;
1634  n_columns++;
1635  }
1636  }
1637  Assert(n_columns <= (int)n_cols, ExcInternalError());
1638 
1639  col_index_ptr = &column_indices[0];
1640  col_value_ptr = &column_values[0];
1641  }
1642 
1643  const int ierr
1644  = MatSetValues (matrix, 1, &petsc_i, n_columns, col_index_ptr,
1645  col_value_ptr, INSERT_VALUES);
1646  AssertThrow (ierr == 0, ExcPETScError(ierr));
1647  }
1648 
1649 
1650 
1651  inline
1652  void
1653  MatrixBase::add (const size_type i,
1654  const size_type j,
1655  const PetscScalar value)
1656  {
1657 
1659 
1660  if (value == PetscScalar())
1661  {
1662  // we have to do checkings on Insert/Add
1663  // in any case
1664  // to be consistent with the MPI
1665  // communication model (see the comments
1666  // in the documentation of
1667  // TrilinosWrappers::Vector), but we can
1668  // save some work if the addend is
1669  // zero. However, these actions are done
1670  // in case we pass on to the other
1671  // function.
1672  prepare_action(LastAction::add);
1673 
1674  return;
1675  }
1676  else
1677  add (i, 1, &j, &value, false);
1678  }
1679 
1680 
1681 
1682  inline
1683  void
1684  MatrixBase::add (const std::vector<size_type> &indices,
1685  const FullMatrix<PetscScalar> &values,
1686  const bool elide_zero_values)
1687  {
1688  Assert (indices.size() == values.m(),
1689  ExcDimensionMismatch(indices.size(), values.m()));
1690  Assert (values.m() == values.n(), ExcNotQuadratic());
1691 
1692  for (size_type i=0; i<indices.size(); ++i)
1693  add (indices[i], indices.size(), &indices[0], &values(i,0),
1694  elide_zero_values);
1695  }
1696 
1697 
1698 
1699  inline
1700  void
1701  MatrixBase::add (const std::vector<size_type> &row_indices,
1702  const std::vector<size_type> &col_indices,
1703  const FullMatrix<PetscScalar> &values,
1704  const bool elide_zero_values)
1705  {
1706  Assert (row_indices.size() == values.m(),
1707  ExcDimensionMismatch(row_indices.size(), values.m()));
1708  Assert (col_indices.size() == values.n(),
1709  ExcDimensionMismatch(col_indices.size(), values.n()));
1710 
1711  for (size_type i=0; i<row_indices.size(); ++i)
1712  add (row_indices[i], col_indices.size(), &col_indices[0], &values(i,0),
1713  elide_zero_values);
1714  }
1715 
1716 
1717 
1718  inline
1719  void
1720  MatrixBase::add (const size_type row,
1721  const std::vector<size_type> &col_indices,
1722  const std::vector<PetscScalar> &values,
1723  const bool elide_zero_values)
1724  {
1725  Assert (col_indices.size() == values.size(),
1726  ExcDimensionMismatch(col_indices.size(), values.size()));
1727 
1728  add (row, col_indices.size(), &col_indices[0], &values[0],
1729  elide_zero_values);
1730  }
1731 
1732 
1733 
1734  inline
1735  void
1736  MatrixBase::add (const size_type row,
1737  const size_type n_cols,
1738  const size_type *col_indices,
1739  const PetscScalar *values,
1740  const bool elide_zero_values,
1741  const bool /*col_indices_are_sorted*/)
1742  {
1743  prepare_action(LastAction::add);
1744 
1745  const PetscInt petsc_i = row;
1746  PetscInt *col_index_ptr;
1747 
1748  PetscScalar const *col_value_ptr;
1749  int n_columns;
1750 
1751  // If we don't elide zeros, the pointers
1752  // are already available...
1753 #ifndef PETSC_USE_64BIT_INDICES
1754  if (elide_zero_values == false)
1755  {
1756  col_index_ptr = (int *)col_indices;
1757  col_value_ptr = values;
1758  n_columns = n_cols;
1759  }
1760  else
1761 #endif
1762  {
1763  // Otherwise, extract nonzero values in
1764  // each row and get the respective index.
1765  if (column_indices.size() < n_cols)
1766  {
1767  column_indices.resize(n_cols);
1768  column_values.resize(n_cols);
1769  }
1770 
1771  n_columns = 0;
1772  for (size_type j=0; j<n_cols; ++j)
1773  {
1774  const PetscScalar value = values[j];
1776  if (value != PetscScalar())
1777  {
1778  column_indices[n_columns] = col_indices[j];
1779  column_values[n_columns] = value;
1780  n_columns++;
1781  }
1782  }
1783  Assert(n_columns <= (int)n_cols, ExcInternalError());
1784 
1785  col_index_ptr = &column_indices[0];
1786  col_value_ptr = &column_values[0];
1787  }
1788 
1789  const int ierr
1790  = MatSetValues (matrix, 1, &petsc_i, n_columns, col_index_ptr,
1791  col_value_ptr, ADD_VALUES);
1792  AssertThrow (ierr == 0, ExcPETScError(ierr));
1793  }
1794 
1795 
1796 
1797 
1798 
1799 
1800  inline
1801  PetscScalar
1802  MatrixBase::operator() (const size_type i,
1803  const size_type j) const
1804  {
1805  return el(i,j);
1806  }
1807 
1808 
1809 
1810  inline
1811  MatrixBase::const_iterator
1812  MatrixBase::begin() const
1813  {
1814  return const_iterator(this, 0, 0);
1815  }
1816 
1817 
1818  inline
1819  MatrixBase::const_iterator
1820  MatrixBase::end() const
1821  {
1822  return const_iterator(this, m(), 0);
1823  }
1824 
1825 
1826  inline
1827  MatrixBase::const_iterator
1828  MatrixBase::begin(const size_type r) const
1829  {
1830  Assert (r < m(), ExcIndexRange(r, 0, m()));
1831  if (row_length(r) > 0)
1832  return const_iterator(this, r, 0);
1833  else
1834  return end (r);
1835  }
1836 
1837 
1838  inline
1839  MatrixBase::const_iterator
1840  MatrixBase::end(const size_type r) const
1841  {
1842  Assert (r < m(), ExcIndexRange(r, 0, m()));
1843 
1844  // place the iterator on the first entry
1845  // past this line, or at the end of the
1846  // matrix
1847  for (size_type i=r+1; i<m(); ++i)
1848  if (row_length(i) > 0)
1849  return const_iterator(this, i, 0);
1850 
1851  // if there is no such line, then take the
1852  // end iterator of the matrix
1853  return end();
1854  }
1855 
1856 
1857 
1858  inline
1859  bool
1860  MatrixBase::in_local_range (const size_type index) const
1861  {
1862  PetscInt begin, end;
1863 
1864  const int ierr = MatGetOwnershipRange (static_cast<const Mat &>(matrix),
1865  &begin, &end);
1866  AssertThrow (ierr == 0, ExcPETScError(ierr));
1867 
1868  return ((index >= static_cast<size_type>(begin)) &&
1869  (index < static_cast<size_type>(end)));
1870  }
1871 
1872 
1873 
1874  inline
1875  void
1876  MatrixBase::prepare_action(const LastAction::Values new_action)
1877  {
1878  if (last_action == new_action)
1879  ;
1880  else if (last_action == LastAction::none)
1881  last_action = new_action;
1882  else
1883  Assert (false, ExcWrongMode (last_action, new_action));
1884  }
1885 
1886 
1887 
1888  inline
1889  void
1890  MatrixBase::prepare_add()
1891  {
1892  prepare_action(LastAction::add);
1893  }
1894 
1895 
1896 
1897  inline
1898  void
1899  MatrixBase::prepare_set()
1900  {
1901  prepare_action(LastAction::insert);
1902  }
1903 
1904 #endif // DOXYGEN
1905 }
1906 
1907 
1908 DEAL_II_NAMESPACE_CLOSE
1909 
1910 
1911 #endif // DEAL_II_WITH_PETSC
1912 
1913 
1914 /*---------------------------- petsc_matrix_base.h ---------------------------*/
1915 
1916 #endif
1917 /*---------------------------- petsc_matrix_base.h ---------------------------*/
std_cxx1x::shared_ptr< const std::vector< size_type > > colnum_cache
PetscScalar matrix_scalar_product(const VectorBase &u, const VectorBase &v) const
void add(const size_type i, const size_type j, const PetscScalar value)
PetscScalar residual(VectorBase &dst, const VectorBase &x, const VectorBase &b) const
size_type n() const
size_type m() const
void prepare_action(const LastAction::Values new_action)
void compress() DEAL_II_DEPRECATED
void Tvmult(VectorBase &dst, const VectorBase &src) const
#define AssertThrow(cond, exc)
Definition: exceptions.h:362
bool operator<(const const_iterator &) const
PetscReal linfty_norm() const
size_type local_size() const
DeclException2(ExcWrongMode, int, int,<< "You tried to do a "<< (arg1==1?"'set'":(arg1==2?"'add'":"???"))<< " operation but the matrix is currently in "<< (arg2==1?"'set'":(arg2==2?"'add'":"???"))<< " mode. You first have to call 'compress()'.")
PetscBool is_symmetric(const double tolerance=1.e-12)
PetscBool is_hermitian(const double tolerance=1.e-12)
const_iterator begin() const
bool is_finite(const double x)
void clear_rows(const std::vector< size_type > &rows, const PetscScalar new_diag_value=0)
std::vector< PetscInt > column_indices
const_iterator(const MatrixBase *matrix, const size_type row, const size_type index)
std::vector< PetscScalar > column_values
MatrixIterators::const_iterator const_iterator
size_type n() const
const_iterator end() const
std::size_t memory_consumption() const
unsigned int global_dof_index
Definition: types.h:100
size_type row_length(const size_type row) const
Accessor(const MatrixBase *matrix, const size_type row, const size_type index)
bool operator==(const const_iterator &) const
#define Assert(cond, exc)
Definition: exceptions.h:299
bool in_local_range(const size_type index) const
bool operator!=(const const_iterator &) const
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
types::global_dof_index size_type
BlockCompressedSparsityPattern CompressedBlockSparsityPattern DEAL_II_DEPRECATED
size_type m() const
PetscReal l1_norm() const
PetscScalar diag_element(const size_type i) const
PetscScalar el(const size_type i, const size_type j) const
::ExceptionBase & ExcIteratorPastEnd()
void write_ascii(const PetscViewerFormat format=PETSC_VIEWER_DEFAULT)
void print(std::ostream &out, const bool alternative_output=false) const
PetscScalar matrix_norm_square(const VectorBase &v) const
DeclException0(ExcSourceEqualsDestination)
::ExceptionBase & ExcNumberNotFinite()
MatrixBase & operator/=(const PetscScalar factor)
std::pair< size_type, size_type > local_range() const
void Tvmult_add(VectorBase &dst, const VectorBase &src) const
LastAction::Values last_action
void set(const size_type i, const size_type j, const PetscScalar value)
DeclException3(ExcAccessToNonlocalRow, int, int, int,<< "You tried to access row "<< arg1<< " of a distributed matrix, but only rows "<< arg2<< " through "<< arg3<< " are stored locally and can be accessed.")
PetscReal frobenius_norm() const
DeclException2(ExcInvalidIndexWithinRow, int, int,<< "Attempt to access element "<< arg2<< " of row "<< arg1<< " which doesn't have that many elements.")
::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
MatrixBase & operator=(const value_type d)
DeclException1(ExcPETScError, int,<< "An error with error number "<< arg1<< " occurred while calling a PETSc function")
void vmult(VectorBase &dst, const VectorBase &src) const
::ExceptionBase & ExcInternalError()
virtual const MPI_Comm & get_mpi_communicator() const =0
std_cxx1x::shared_ptr< const std::vector< PetscScalar > > value_cache
size_type n_nonzero_elements() const
void clear_row(const size_type row, const PetscScalar new_diag_value=0)
MatrixBase & operator*=(const PetscScalar factor)
void vmult_add(VectorBase &dst, const VectorBase &src) const