Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
petsc_matrix_base.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_petsc_matrix_base_h
17 # define dealii_petsc_matrix_base_h
18 
19 
20 # include <deal.II/base/config.h>
21 
22 # ifdef DEAL_II_WITH_PETSC
23 
24 # include <deal.II/base/subscriptor.h>
25 
26 # include <deal.II/lac/exceptions.h>
27 # include <deal.II/lac/full_matrix.h>
28 # include <deal.II/lac/petsc_compatibility.h>
29 # include <deal.II/lac/petsc_vector_base.h>
30 # include <deal.II/lac/vector_operation.h>
31 
32 # include <petscmat.h>
33 
34 # include <cmath>
35 # include <memory>
36 # include <vector>
37 
38 DEAL_II_NAMESPACE_OPEN
39 
40 template <typename Matrix>
41 class BlockMatrixBase;
42 
43 
44 namespace PETScWrappers
45 {
46  // forward declarations
47  class MatrixBase;
48 
49  namespace MatrixIterators
50  {
66  {
67  private:
71  class Accessor
72  {
73  public:
78 
83  Accessor(const MatrixBase *matrix,
84  const size_type row,
85  const size_type index);
86 
90  size_type
91  row() const;
92 
96  size_type
97  index() const;
98 
102  size_type
103  column() const;
104 
108  PetscScalar
109  value() const;
110 
119  int,
120  int,
121  int,
122  << "You tried to access row " << arg1
123  << " of a distributed matrix, but only rows " << arg2
124  << " through " << arg3
125  << " are stored locally and can be accessed.");
126 
127  private:
131  mutable MatrixBase *matrix;
132 
137 
142 
155  std::shared_ptr<const std::vector<size_type>> colnum_cache;
156 
160  std::shared_ptr<const std::vector<PetscScalar>> value_cache;
161 
167  void
169 
173  friend class const_iterator;
174  };
175 
176  public:
181 
186  const_iterator(const MatrixBase *matrix,
187  const size_type row,
188  const size_type index);
189 
194  operator++();
195 
200  operator++(int);
201 
205  const Accessor &operator*() const;
206 
210  const Accessor *operator->() const;
211 
216  bool
217  operator==(const const_iterator &) const;
221  bool
222  operator!=(const const_iterator &) const;
223 
229  bool
230  operator<(const const_iterator &) const;
231 
236  int,
237  int,
238  << "Attempt to access element " << arg2 << " of row "
239  << arg1 << " which doesn't have that many elements.");
240 
241  private:
246  };
247 
248  } // namespace MatrixIterators
249 
250 
283  class MatrixBase : public Subscriptor
284  {
285  public:
290 
295 
299  using value_type = PetscScalar;
300 
304  MatrixBase();
305 
311  MatrixBase(const MatrixBase &) = delete;
312 
318  MatrixBase &
319  operator=(const MatrixBase &) = delete;
320 
324  virtual ~MatrixBase() override;
325 
335  MatrixBase &
336  operator=(const value_type d);
341  void
342  clear();
343 
353  void
354  set(const size_type i, const size_type j, const PetscScalar value);
355 
375  void
376  set(const std::vector<size_type> & indices,
377  const FullMatrix<PetscScalar> &full_matrix,
378  const bool elide_zero_values = false);
379 
385  void
386  set(const std::vector<size_type> & row_indices,
387  const std::vector<size_type> & col_indices,
388  const FullMatrix<PetscScalar> &full_matrix,
389  const bool elide_zero_values = false);
390 
405  void
406  set(const size_type row,
407  const std::vector<size_type> & col_indices,
408  const std::vector<PetscScalar> &values,
409  const bool elide_zero_values = false);
410 
425  void
426  set(const size_type row,
427  const size_type n_cols,
428  const size_type * col_indices,
429  const PetscScalar *values,
430  const bool elide_zero_values = false);
431 
441  void
442  add(const size_type i, const size_type j, const PetscScalar value);
443 
463  void
464  add(const std::vector<size_type> & indices,
465  const FullMatrix<PetscScalar> &full_matrix,
466  const bool elide_zero_values = true);
467 
473  void
474  add(const std::vector<size_type> & row_indices,
475  const std::vector<size_type> & col_indices,
476  const FullMatrix<PetscScalar> &full_matrix,
477  const bool elide_zero_values = true);
478 
493  void
494  add(const size_type row,
495  const std::vector<size_type> & col_indices,
496  const std::vector<PetscScalar> &values,
497  const bool elide_zero_values = true);
498 
513  void
514  add(const size_type row,
515  const size_type n_cols,
516  const size_type * col_indices,
517  const PetscScalar *values,
518  const bool elide_zero_values = true,
519  const bool col_indices_are_sorted = false);
520 
537  void
538  clear_row(const size_type row, const PetscScalar new_diag_value = 0);
539 
548  void
549  clear_rows(const std::vector<size_type> &rows,
550  const PetscScalar new_diag_value = 0);
551 
563  void
564  compress(const VectorOperation::values operation);
565 
577  PetscScalar
578  operator()(const size_type i, const size_type j) const;
579 
587  PetscScalar
588  el(const size_type i, const size_type j) const;
589 
599  PetscScalar
600  diag_element(const size_type i) const;
601 
605  size_type
606  m() const;
607 
611  size_type
612  n() const;
613 
622  size_type
623  local_size() const;
624 
633  std::pair<size_type, size_type>
634  local_range() const;
635 
640  bool
641  in_local_range(const size_type index) const;
642 
647  virtual const MPI_Comm &
648  get_mpi_communicator() const = 0;
649 
655  size_type
656  n_nonzero_elements() const;
657 
661  size_type
662  row_length(const size_type row) const;
663 
671  PetscReal
672  l1_norm() const;
673 
681  PetscReal
682  linfty_norm() const;
683 
688  PetscReal
689  frobenius_norm() const;
690 
691 
711  PetscScalar
712  matrix_norm_square(const VectorBase &v) const;
713 
714 
728  PetscScalar
729  matrix_scalar_product(const VectorBase &u, const VectorBase &v) const;
730 
735  PetscScalar
736  trace() const;
737 
741  MatrixBase &
742  operator*=(const PetscScalar factor);
743 
747  MatrixBase &
748  operator/=(const PetscScalar factor);
749 
750 
755  MatrixBase &
756  add(const PetscScalar factor, const MatrixBase &other);
757 
758 
764  DEAL_II_DEPRECATED
765  MatrixBase &
766  add(const MatrixBase &other, const PetscScalar factor);
767 
779  void
780  vmult(VectorBase &dst, const VectorBase &src) const;
781 
794  void
795  Tvmult(VectorBase &dst, const VectorBase &src) const;
796 
808  void
809  vmult_add(VectorBase &dst, const VectorBase &src) const;
810 
823  void
824  Tvmult_add(VectorBase &dst, const VectorBase &src) const;
825 
838  PetscScalar
839  residual(VectorBase &dst, const VectorBase &x, const VectorBase &b) const;
840 
845  begin() const;
846 
851  end() const;
852 
862  begin(const size_type r) const;
863 
873  end(const size_type r) const;
874 
882  operator Mat() const;
883 
889  Mat &
890  petsc_matrix();
891 
895  void
896  transpose();
897 
902  PetscBool
903  is_symmetric(const double tolerance = 1.e-12);
904 
910  PetscBool
911  is_hermitian(const double tolerance = 1.e-12);
912 
919  void
920  write_ascii(const PetscViewerFormat format = PETSC_VIEWER_DEFAULT);
921 
929  void
930  print(std::ostream &out, const bool alternative_output = false) const;
931 
935  std::size_t
936  memory_consumption() const;
937 
942  "You are attempting an operation on two matrices that "
943  "are the same object, but the operation requires that the "
944  "two objects are in fact different.");
945 
950  int,
951  int,
952  << "You tried to do a "
953  << (arg1 == 1 ? "'set'" : (arg1 == 2 ? "'add'" : "???"))
954  << " operation but the matrix is currently in "
955  << (arg2 == 1 ? "'set'" : (arg2 == 2 ? "'add'" : "???"))
956  << " mode. You first have to call 'compress()'.");
957 
958  protected:
963  Mat matrix;
964 
969 
975  void
976  prepare_action(const VectorOperation::values new_action);
977 
983  void
985 
996  void
997  prepare_add();
1003  void
1004  prepare_set();
1005 
1021  void
1022  mmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const;
1023 
1040  void
1041  Tmmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const;
1042 
1043  private:
1058  mutable std::vector<PetscInt> column_indices;
1059 
1068  mutable std::vector<PetscScalar> column_values;
1069 
1070 
1074  template <class>
1075  friend class ::BlockMatrixBase;
1076  };
1077 
1078 
1079 
1080 # ifndef DOXYGEN
1081  // ---------------------- inline and template functions ---------------------
1082 
1083 
1084  namespace MatrixIterators
1085  {
1086  inline const_iterator::Accessor::Accessor(const MatrixBase *matrix,
1087  const size_type row,
1088  const size_type index)
1089  : matrix(const_cast<MatrixBase *>(matrix))
1090  , a_row(row)
1091  , a_index(index)
1092  {
1093  visit_present_row();
1094  }
1095 
1096 
1097 
1100  {
1101  Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1102  return a_row;
1103  }
1104 
1105 
1108  {
1109  Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1110  return (*colnum_cache)[a_index];
1111  }
1112 
1113 
1116  {
1117  Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1118  return a_index;
1119  }
1120 
1121 
1122  inline PetscScalar
1124  {
1125  Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1126  return (*value_cache)[a_index];
1127  }
1128 
1129 
1130  inline const_iterator::const_iterator(const MatrixBase *matrix,
1131  const size_type row,
1132  const size_type index)
1133  : accessor(matrix, row, index)
1134  {}
1135 
1136 
1137 
1138  inline const_iterator &
1139  const_iterator::operator++()
1140  {
1141  Assert(accessor.a_row < accessor.matrix->m(), ExcIteratorPastEnd());
1142 
1143  ++accessor.a_index;
1144 
1145  // if at end of line: do one step, then cycle until we find a
1146  // row with a nonzero number of entries
1147  if (accessor.a_index >= accessor.colnum_cache->size())
1148  {
1149  accessor.a_index = 0;
1150  ++accessor.a_row;
1151 
1152  while ((accessor.a_row < accessor.matrix->m()) &&
1153  (accessor.a_row < accessor.matrix->local_range().second) &&
1154  (accessor.matrix->row_length(accessor.a_row) == 0))
1155  ++accessor.a_row;
1156 
1157  accessor.visit_present_row();
1158  }
1159  return *this;
1160  }
1161 
1162 
1163  inline const_iterator
1164  const_iterator::operator++(int)
1165  {
1166  const const_iterator old_state = *this;
1167  ++(*this);
1168  return old_state;
1169  }
1170 
1171 
1172  inline const const_iterator::Accessor &const_iterator::operator*() const
1173  {
1174  return accessor;
1175  }
1176 
1177 
1178  inline const const_iterator::Accessor *const_iterator::operator->() const
1179  {
1180  return &accessor;
1181  }
1182 
1183 
1184  inline bool
1185  const_iterator::operator==(const const_iterator &other) const
1186  {
1187  return (accessor.a_row == other.accessor.a_row &&
1188  accessor.a_index == other.accessor.a_index);
1189  }
1190 
1191 
1192  inline bool
1193  const_iterator::operator!=(const const_iterator &other) const
1194  {
1195  return !(*this == other);
1196  }
1197 
1198 
1199  inline bool
1200  const_iterator::operator<(const const_iterator &other) const
1201  {
1202  return (accessor.row() < other.accessor.row() ||
1203  (accessor.row() == other.accessor.row() &&
1204  accessor.index() < other.accessor.index()));
1205  }
1206 
1207  } // namespace MatrixIterators
1208 
1209 
1210 
1211  // Inline the set() and add()
1212  // functions, since they will be
1213  // called frequently, and the
1214  // compiler can optimize away
1215  // some unnecessary loops when
1216  // the sizes are given at
1217  // compile time.
1218  inline void
1219  MatrixBase::set(const size_type i, const size_type j, const PetscScalar value)
1220  {
1221  AssertIsFinite(value);
1222 
1223  set(i, 1, &j, &value, false);
1224  }
1225 
1226 
1227 
1228  inline void
1229  MatrixBase::set(const std::vector<size_type> & indices,
1230  const FullMatrix<PetscScalar> &values,
1231  const bool elide_zero_values)
1232  {
1233  Assert(indices.size() == values.m(),
1234  ExcDimensionMismatch(indices.size(), values.m()));
1235  Assert(values.m() == values.n(), ExcNotQuadratic());
1236 
1237  for (size_type i = 0; i < indices.size(); ++i)
1238  set(indices[i],
1239  indices.size(),
1240  indices.data(),
1241  &values(i, 0),
1242  elide_zero_values);
1243  }
1244 
1245 
1246 
1247  inline void
1248  MatrixBase::set(const std::vector<size_type> & row_indices,
1249  const std::vector<size_type> & col_indices,
1250  const FullMatrix<PetscScalar> &values,
1251  const bool elide_zero_values)
1252  {
1253  Assert(row_indices.size() == values.m(),
1254  ExcDimensionMismatch(row_indices.size(), values.m()));
1255  Assert(col_indices.size() == values.n(),
1256  ExcDimensionMismatch(col_indices.size(), values.n()));
1257 
1258  for (size_type i = 0; i < row_indices.size(); ++i)
1259  set(row_indices[i],
1260  col_indices.size(),
1261  col_indices.data(),
1262  &values(i, 0),
1263  elide_zero_values);
1264  }
1265 
1266 
1267 
1268  inline void
1269  MatrixBase::set(const size_type row,
1270  const std::vector<size_type> & col_indices,
1271  const std::vector<PetscScalar> &values,
1272  const bool elide_zero_values)
1273  {
1274  Assert(col_indices.size() == values.size(),
1275  ExcDimensionMismatch(col_indices.size(), values.size()));
1276 
1277  set(row,
1278  col_indices.size(),
1279  col_indices.data(),
1280  values.data(),
1281  elide_zero_values);
1282  }
1283 
1284 
1285 
1286  inline void
1287  MatrixBase::set(const size_type row,
1288  const size_type n_cols,
1289  const size_type * col_indices,
1290  const PetscScalar *values,
1291  const bool elide_zero_values)
1292  {
1293  prepare_action(VectorOperation::insert);
1294 
1295  const PetscInt petsc_i = row;
1296  PetscInt const *col_index_ptr;
1297 
1298  PetscScalar const *col_value_ptr;
1299  int n_columns;
1300 
1301  // If we don't elide zeros, the pointers are already available...
1302  if (elide_zero_values == false)
1303  {
1304  col_index_ptr = reinterpret_cast<const PetscInt *>(col_indices);
1305  col_value_ptr = values;
1306  n_columns = n_cols;
1307  }
1308  else
1309  {
1310  // Otherwise, extract nonzero values in each row and get the
1311  // respective index.
1312  if (column_indices.size() < n_cols)
1313  {
1314  column_indices.resize(n_cols);
1315  column_values.resize(n_cols);
1316  }
1317 
1318  n_columns = 0;
1319  for (size_type j = 0; j < n_cols; ++j)
1320  {
1321  const PetscScalar value = values[j];
1322  AssertIsFinite(value);
1323  if (value != PetscScalar())
1324  {
1325  column_indices[n_columns] = col_indices[j];
1326  column_values[n_columns] = value;
1327  n_columns++;
1328  }
1329  }
1330  AssertIndexRange(n_columns, n_cols + 1);
1331 
1332  col_index_ptr = column_indices.data();
1333  col_value_ptr = column_values.data();
1334  }
1335 
1336  const PetscErrorCode ierr = MatSetValues(matrix,
1337  1,
1338  &petsc_i,
1339  n_columns,
1340  col_index_ptr,
1341  col_value_ptr,
1342  INSERT_VALUES);
1343  AssertThrow(ierr == 0, ExcPETScError(ierr));
1344  }
1345 
1346 
1347 
1348  inline void
1349  MatrixBase::add(const size_type i, const size_type j, const PetscScalar value)
1350  {
1351  AssertIsFinite(value);
1352 
1353  if (value == PetscScalar())
1354  {
1355  // we have to check after using Insert/Add in any case to be
1356  // consistent with the MPI communication model, but we can save
1357  // some work if the addend is zero. However, these actions are done
1358  // in case we pass on to the other function.
1359  prepare_action(VectorOperation::add);
1360 
1361  return;
1362  }
1363  else
1364  add(i, 1, &j, &value, false);
1365  }
1366 
1367 
1368 
1369  inline void
1370  MatrixBase::add(const std::vector<size_type> & indices,
1371  const FullMatrix<PetscScalar> &values,
1372  const bool elide_zero_values)
1373  {
1374  Assert(indices.size() == values.m(),
1375  ExcDimensionMismatch(indices.size(), values.m()));
1376  Assert(values.m() == values.n(), ExcNotQuadratic());
1377 
1378  for (size_type i = 0; i < indices.size(); ++i)
1379  add(indices[i],
1380  indices.size(),
1381  indices.data(),
1382  &values(i, 0),
1383  elide_zero_values);
1384  }
1385 
1386 
1387 
1388  inline void
1389  MatrixBase::add(const std::vector<size_type> & row_indices,
1390  const std::vector<size_type> & col_indices,
1391  const FullMatrix<PetscScalar> &values,
1392  const bool elide_zero_values)
1393  {
1394  Assert(row_indices.size() == values.m(),
1395  ExcDimensionMismatch(row_indices.size(), values.m()));
1396  Assert(col_indices.size() == values.n(),
1397  ExcDimensionMismatch(col_indices.size(), values.n()));
1398 
1399  for (size_type i = 0; i < row_indices.size(); ++i)
1400  add(row_indices[i],
1401  col_indices.size(),
1402  col_indices.data(),
1403  &values(i, 0),
1404  elide_zero_values);
1405  }
1406 
1407 
1408 
1409  inline void
1410  MatrixBase::add(const size_type row,
1411  const std::vector<size_type> & col_indices,
1412  const std::vector<PetscScalar> &values,
1413  const bool elide_zero_values)
1414  {
1415  Assert(col_indices.size() == values.size(),
1416  ExcDimensionMismatch(col_indices.size(), values.size()));
1417 
1418  add(row,
1419  col_indices.size(),
1420  col_indices.data(),
1421  values.data(),
1422  elide_zero_values);
1423  }
1424 
1425 
1426 
1427  inline void
1428  MatrixBase::add(const size_type row,
1429  const size_type n_cols,
1430  const size_type * col_indices,
1431  const PetscScalar *values,
1432  const bool elide_zero_values,
1433  const bool /*col_indices_are_sorted*/)
1434  {
1435  (void)elide_zero_values;
1436 
1437  prepare_action(VectorOperation::add);
1438 
1439  const PetscInt petsc_i = row;
1440  PetscInt const *col_index_ptr;
1441 
1442  PetscScalar const *col_value_ptr;
1443  int n_columns;
1444 
1445  // If we don't elide zeros, the pointers are already available...
1446  if (elide_zero_values == false)
1447  {
1448  col_index_ptr = reinterpret_cast<const PetscInt *>(col_indices);
1449  col_value_ptr = values;
1450  n_columns = n_cols;
1451  }
1452  else
1453  {
1454  // Otherwise, extract nonzero values in each row and get the
1455  // respective index.
1456  if (column_indices.size() < n_cols)
1457  {
1458  column_indices.resize(n_cols);
1459  column_values.resize(n_cols);
1460  }
1461 
1462  n_columns = 0;
1463  for (size_type j = 0; j < n_cols; ++j)
1464  {
1465  const PetscScalar value = values[j];
1466  AssertIsFinite(value);
1467  if (value != PetscScalar())
1468  {
1469  column_indices[n_columns] = col_indices[j];
1470  column_values[n_columns] = value;
1471  n_columns++;
1472  }
1473  }
1474  AssertIndexRange(n_columns, n_cols + 1);
1475 
1476  col_index_ptr = column_indices.data();
1477  col_value_ptr = column_values.data();
1478  }
1479 
1480  const PetscErrorCode ierr = MatSetValues(
1481  matrix, 1, &petsc_i, n_columns, col_index_ptr, col_value_ptr, ADD_VALUES);
1482  AssertThrow(ierr == 0, ExcPETScError(ierr));
1483  }
1484 
1485 
1486 
1487  inline PetscScalar
1488  MatrixBase::operator()(const size_type i, const size_type j) const
1489  {
1490  return el(i, j);
1491  }
1492 
1493 
1494 
1495  inline MatrixBase::const_iterator
1496  MatrixBase::begin() const
1497  {
1498  return const_iterator(this, 0, 0);
1499  }
1500 
1501 
1502  inline MatrixBase::const_iterator
1503  MatrixBase::end() const
1504  {
1505  return const_iterator(this, m(), 0);
1506  }
1507 
1508 
1509  inline MatrixBase::const_iterator
1510  MatrixBase::begin(const size_type r) const
1511  {
1512  Assert(in_local_range(r),
1513  ExcIndexRange(r, local_range().first, local_range().second));
1514 
1515  if (row_length(r) > 0)
1516  return const_iterator(this, r, 0);
1517  else
1518  return end(r);
1519  }
1520 
1521 
1522  inline MatrixBase::const_iterator
1523  MatrixBase::end(const size_type r) const
1524  {
1525  Assert(in_local_range(r),
1526  ExcIndexRange(r, local_range().first, local_range().second));
1527 
1528  // place the iterator on the first entry past this line, or at the
1529  // end of the matrix
1530  //
1531  // in the parallel case, we need to put it on the first entry of
1532  // the first row after the locally owned range. this of course
1533  // doesn't exist, but we can nevertheless create such an
1534  // iterator. we need to check whether 'i' is past the locally
1535  // owned range of rows first, before we ask for the length of the
1536  // row since the latter query leads to an exception in case we ask
1537  // for a row that is not locally owned
1538  for (size_type i = r + 1; i < m(); ++i)
1539  if (i == local_range().second || (row_length(i) > 0))
1540  return const_iterator(this, i, 0);
1541 
1542  // if there is no such line, then take the
1543  // end iterator of the matrix
1544  return end();
1545  }
1546 
1547 
1548 
1549  inline bool
1550  MatrixBase::in_local_range(const size_type index) const
1551  {
1552  PetscInt begin, end;
1553 
1554  const PetscErrorCode ierr =
1555  MatGetOwnershipRange(static_cast<const Mat &>(matrix), &begin, &end);
1556  AssertThrow(ierr == 0, ExcPETScError(ierr));
1557 
1558  return ((index >= static_cast<size_type>(begin)) &&
1559  (index < static_cast<size_type>(end)));
1560  }
1561 
1562 
1563 
1564  inline void
1565  MatrixBase::prepare_action(const VectorOperation::values new_action)
1566  {
1567  if (last_action == VectorOperation::unknown)
1568  last_action = new_action;
1569 
1570  Assert(last_action == new_action, ExcWrongMode(last_action, new_action));
1571  }
1572 
1573 
1574 
1575  inline void
1576  MatrixBase::assert_is_compressed()
1577  {
1578  // compress() sets the last action to none, which allows us to check if
1579  // there are pending add/insert operations:
1580  AssertThrow(last_action == VectorOperation::unknown,
1581  ExcMessage("Error: missing compress() call."));
1582  }
1583 
1584 
1585 
1586  inline void
1587  MatrixBase::prepare_add()
1588  {
1589  prepare_action(VectorOperation::add);
1590  }
1591 
1592 
1593 
1594  inline void
1595  MatrixBase::prepare_set()
1596  {
1597  prepare_action(VectorOperation::insert);
1598  }
1599 
1600 # endif // DOXYGEN
1601 } // namespace PETScWrappers
1602 
1603 
1604 DEAL_II_NAMESPACE_CLOSE
1605 
1606 
1607 # endif // DEAL_II_WITH_PETSC
1608 
1609 #endif
1610 /*---------------------------- petsc_matrix_base.h --------------------------*/
size_type m() const
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:541
bool operator==(const const_iterator &) const
void Tmmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const
Contents is actually a matrix.
PetscScalar matrix_scalar_product(const VectorBase &u, const VectorBase &v) const
static ::ExceptionBase & ExcBeyondEndOfMatrix()
void add(const size_type i, const size_type j, const PetscScalar value)
void vmult(VectorBase &dst, const VectorBase &src) const
const_iterator begin() const
PetscReal frobenius_norm() const
static ::ExceptionBase & ExcInvalidIndexWithinRow(int arg1, int arg2)
std::shared_ptr< const std::vector< PetscScalar > > value_cache
#define AssertIndexRange(index, range)
Definition: exceptions.h:1637
size_type row_length(const size_type row) const
#define AssertThrow(cond, exc)
Definition: exceptions.h:1519
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
SymmetricTensor< rank_, dim, Number > operator*(const SymmetricTensor< rank_, dim, Number > &t, const Number &factor)
PetscBool is_symmetric(const double tolerance=1.e-12)
PetscBool is_hermitian(const double tolerance=1.e-12)
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)
void compress(const VectorOperation::values operation)
std::vector< PetscScalar > column_values
size_type n_nonzero_elements() const
MatrixBase & operator=(const MatrixBase &)=delete
static ::ExceptionBase & ExcMessage(std::string arg1)
size_type n() 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:1407
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
#define DeclException0(Exception0)
Definition: exceptions.h:473
void vmult_add(VectorBase &dst, const VectorBase &src) const
PetscScalar residual(VectorBase &dst, const VectorBase &x, const VectorBase &b) const
void mmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const
PetscScalar diag_element(const size_type i) const
void Tvmult(VectorBase &dst, const VectorBase &src) const
PetscScalar el(const size_type i, const size_type j) const
void Tvmult_add(VectorBase &dst, const VectorBase &src) const
const_iterator end() const
PetscScalar operator()(const size_type i, const size_type j) const
std::shared_ptr< const std::vector< size_type > > colnum_cache
static ::ExceptionBase & ExcIteratorPastEnd()
void print(std::ostream &out, const bool alternative_output=false) const
static ::ExceptionBase & ExcNotQuadratic()
virtual ~MatrixBase() override
void write_ascii(const PetscViewerFormat format=PETSC_VIEWER_DEFAULT)
VectorOperation::values last_action
void prepare_action(const VectorOperation::values new_action)
unsigned int global_dof_index
Definition: types.h:89
__global__ void set(Number *val, const Number s, const size_type N)
types::global_dof_index size_type
bool in_local_range(const size_type index) const
MatrixBase & operator*=(const PetscScalar factor)
std::pair< size_type, size_type > local_range() const
static ::ExceptionBase & ExcAccessToNonlocalRow(int arg1, int arg2, int arg3)
PetscScalar matrix_norm_square(const VectorBase &v) const
MatrixBase & operator/=(const PetscScalar factor)
static ::ExceptionBase & ExcSourceEqualsDestination()
static ::ExceptionBase & ExcWrongMode(int arg1, int arg2)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition: exceptions.h:564
bool operator!=(const const_iterator &) const
virtual const MPI_Comm & get_mpi_communicator() const =0
#define AssertIsFinite(number)
Definition: exceptions.h:1669
void clear_row(const size_type row, const PetscScalar new_diag_value=0)
std::size_t memory_consumption() const