Reference documentation for deal.II version Git 9e557027ad 2021-09-25 18:07:42 +0200
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
petsc_matrix_base.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2020 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 
25 
26 # include <deal.II/lac/exceptions.h>
27 # include <deal.II/lac/full_matrix.h>
31 
32 # include <petscmat.h>
33 
34 # include <cmath>
35 # include <memory>
36 # include <vector>
37 
39 
40 // Forward declarations
41 # ifndef DOXYGEN
42 template <typename Matrix>
43 class BlockMatrixBase;
44 # endif
45 
46 
47 namespace PETScWrappers
48 {
49  // forward declarations
50  class MatrixBase;
51 
52  namespace MatrixIterators
53  {
68  {
69  private:
73  class Accessor
74  {
75  public:
80 
85  Accessor(const MatrixBase *matrix,
86  const size_type row,
87  const size_type index);
88 
92  size_type
93  row() const;
94 
98  size_type
99  index() const;
100 
104  size_type
105  column() const;
106 
110  PetscScalar
111  value() const;
112 
121  int,
122  int,
123  int,
124  << "You tried to access row " << arg1
125  << " of a distributed matrix, but only rows " << arg2
126  << " through " << arg3
127  << " are stored locally and can be accessed.");
128 
129  private:
133  mutable MatrixBase *matrix;
134 
139 
144 
157  std::shared_ptr<const std::vector<size_type>> colnum_cache;
158 
162  std::shared_ptr<const std::vector<PetscScalar>> value_cache;
163 
169  void
171 
172  // Make enclosing class a friend.
173  friend class const_iterator;
174  };
175 
176  public:
181 
187  const size_type row,
188  const size_type index);
189 
194  operator++();
195 
200  operator++(int);
201 
205  const Accessor &
206  operator*() const;
207 
211  const Accessor *
212  operator->() const;
213 
218  bool
219  operator==(const const_iterator &) const;
223  bool
224  operator!=(const const_iterator &) const;
225 
231  bool
232  operator<(const const_iterator &) const;
233 
238  int,
239  int,
240  << "Attempt to access element " << arg2 << " of row "
241  << arg1 << " which doesn't have that many elements.");
242 
243  private:
248  };
249 
250  } // namespace MatrixIterators
251 
252 
284  class MatrixBase : public Subscriptor
285  {
286  public:
291 
296 
300  using value_type = PetscScalar;
301 
305  MatrixBase();
306 
312  MatrixBase(const MatrixBase &) = delete;
313 
319  MatrixBase &
320  operator=(const MatrixBase &) = delete;
321 
325  virtual ~MatrixBase() override;
326 
336  MatrixBase &
337  operator=(const value_type d);
342  void
343  clear();
344 
354  void
355  set(const size_type i, const size_type j, const PetscScalar value);
356 
376  void
377  set(const std::vector<size_type> & indices,
378  const FullMatrix<PetscScalar> &full_matrix,
379  const bool elide_zero_values = false);
380 
386  void
387  set(const std::vector<size_type> & row_indices,
388  const std::vector<size_type> & col_indices,
389  const FullMatrix<PetscScalar> &full_matrix,
390  const bool elide_zero_values = false);
391 
406  void
407  set(const size_type row,
408  const std::vector<size_type> & col_indices,
409  const std::vector<PetscScalar> &values,
410  const bool elide_zero_values = false);
411 
426  void
427  set(const size_type row,
428  const size_type n_cols,
429  const size_type * col_indices,
430  const PetscScalar *values,
431  const bool elide_zero_values = false);
432 
442  void
443  add(const size_type i, const size_type j, const PetscScalar value);
444 
464  void
465  add(const std::vector<size_type> & indices,
466  const FullMatrix<PetscScalar> &full_matrix,
467  const bool elide_zero_values = true);
468 
474  void
475  add(const std::vector<size_type> & row_indices,
476  const std::vector<size_type> & col_indices,
477  const FullMatrix<PetscScalar> &full_matrix,
478  const bool elide_zero_values = true);
479 
494  void
495  add(const size_type row,
496  const std::vector<size_type> & col_indices,
497  const std::vector<PetscScalar> &values,
498  const bool elide_zero_values = true);
499 
514  void
515  add(const size_type row,
516  const size_type n_cols,
517  const size_type * col_indices,
518  const PetscScalar *values,
519  const bool elide_zero_values = true,
520  const bool col_indices_are_sorted = false);
521 
538  void
539  clear_row(const size_type row, const PetscScalar new_diag_value = 0);
540 
549  void
550  clear_rows(const std::vector<size_type> &rows,
551  const PetscScalar new_diag_value = 0);
552 
564  void
565  compress(const VectorOperation::values operation);
566 
578  PetscScalar
579  operator()(const size_type i, const size_type j) const;
580 
588  PetscScalar
589  el(const size_type i, const size_type j) const;
590 
600  PetscScalar
601  diag_element(const size_type i) const;
602 
606  size_type
607  m() const;
608 
612  size_type
613  n() const;
614 
623  size_type
624  local_size() const;
625 
634  std::pair<size_type, size_type>
635  local_range() const;
636 
641  bool
642  in_local_range(const size_type index) const;
643 
648  virtual const MPI_Comm &
649  get_mpi_communicator() const = 0;
650 
656  size_type
657  n_nonzero_elements() const;
658 
662  size_type
663  row_length(const size_type row) const;
664 
672  PetscReal
673  l1_norm() const;
674 
682  PetscReal
683  linfty_norm() const;
684 
689  PetscReal
690  frobenius_norm() const;
691 
692 
712  PetscScalar
713  matrix_norm_square(const VectorBase &v) const;
714 
715 
729  PetscScalar
730  matrix_scalar_product(const VectorBase &u, const VectorBase &v) const;
731 
736  PetscScalar
737  trace() const;
738 
742  MatrixBase &
743  operator*=(const PetscScalar factor);
744 
748  MatrixBase &
749  operator/=(const PetscScalar factor);
750 
751 
756  MatrixBase &
757  add(const PetscScalar factor, const MatrixBase &other);
758 
770  void
771  vmult(VectorBase &dst, const VectorBase &src) const;
772 
785  void
786  Tvmult(VectorBase &dst, const VectorBase &src) const;
787 
799  void
800  vmult_add(VectorBase &dst, const VectorBase &src) const;
801 
814  void
815  Tvmult_add(VectorBase &dst, const VectorBase &src) const;
816 
829  PetscScalar
830  residual(VectorBase &dst, const VectorBase &x, const VectorBase &b) const;
831 
838  begin() const;
839 
846  end() const;
847 
857  begin(const size_type r) const;
858 
868  end(const size_type r) const;
869 
877  operator Mat() const;
878 
884  Mat &
885  petsc_matrix();
886 
890  void
891  transpose();
892 
897  PetscBool
898  is_symmetric(const double tolerance = 1.e-12);
899 
905  PetscBool
906  is_hermitian(const double tolerance = 1.e-12);
907 
914  void
915  write_ascii(const PetscViewerFormat format = PETSC_VIEWER_DEFAULT);
916 
924  void
925  print(std::ostream &out, const bool alternative_output = false) const;
926 
930  std::size_t
931  memory_consumption() const;
932 
936  DeclExceptionMsg(ExcSourceEqualsDestination,
937  "You are attempting an operation on two matrices that "
938  "are the same object, but the operation requires that the "
939  "two objects are in fact different.");
940 
944  DeclException2(ExcWrongMode,
945  int,
946  int,
947  << "You tried to do a "
948  << (arg1 == 1 ? "'set'" : (arg1 == 2 ? "'add'" : "???"))
949  << " operation but the matrix is currently in "
950  << (arg2 == 1 ? "'set'" : (arg2 == 2 ? "'add'" : "???"))
951  << " mode. You first have to call 'compress()'.");
952 
953  protected:
958  Mat matrix;
959 
964 
970  void
971  prepare_action(const VectorOperation::values new_action);
972 
978  void
979  assert_is_compressed();
980 
991  void
992  prepare_add();
998  void
999  prepare_set();
1000 
1016  void
1017  mmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const;
1018 
1035  void
1036  Tmmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const;
1037 
1038  private:
1053  mutable std::vector<PetscInt> column_indices;
1054 
1063  mutable std::vector<PetscScalar> column_values;
1064 
1065 
1066  // To allow calling protected prepare_add() and prepare_set().
1067  template <class>
1068  friend class ::BlockMatrixBase;
1069  };
1070 
1071 
1072 
1073 # ifndef DOXYGEN
1074  // ---------------------- inline and template functions ---------------------
1075 
1076 
1077  namespace MatrixIterators
1078  {
1080  const size_type row,
1081  const size_type index)
1082  : matrix(const_cast<MatrixBase *>(matrix))
1083  , a_row(row)
1084  , a_index(index)
1085  {
1087  }
1088 
1089 
1090 
1093  {
1094  Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1095  return a_row;
1096  }
1097 
1098 
1101  {
1102  Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1103  return (*colnum_cache)[a_index];
1104  }
1105 
1106 
1109  {
1110  Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1111  return a_index;
1112  }
1113 
1114 
1115  inline PetscScalar
1117  {
1118  Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1119  return (*value_cache)[a_index];
1120  }
1121 
1122 
1123  inline const_iterator::const_iterator(const MatrixBase *matrix,
1124  const size_type row,
1125  const size_type index)
1126  : accessor(matrix, row, index)
1127  {}
1128 
1129 
1130 
1131  inline const_iterator &
1133  {
1135 
1136  ++accessor.a_index;
1137 
1138  // if at end of line: do one step, then cycle until we find a
1139  // row with a nonzero number of entries
1140  if (accessor.a_index >= accessor.colnum_cache->size())
1141  {
1142  accessor.a_index = 0;
1143  ++accessor.a_row;
1144 
1145  while ((accessor.a_row < accessor.matrix->m()) &&
1146  (accessor.a_row < accessor.matrix->local_range().second) &&
1148  ++accessor.a_row;
1149 
1151  }
1152  return *this;
1153  }
1154 
1155 
1156  inline const_iterator
1158  {
1159  const const_iterator old_state = *this;
1160  ++(*this);
1161  return old_state;
1162  }
1163 
1164 
1165  inline const const_iterator::Accessor &
1167  {
1168  return accessor;
1169  }
1170 
1171 
1172  inline const const_iterator::Accessor *
1174  {
1175  return &accessor;
1176  }
1177 
1178 
1179  inline bool
1180  const_iterator::operator==(const const_iterator &other) const
1181  {
1182  return (accessor.a_row == other.accessor.a_row &&
1183  accessor.a_index == other.accessor.a_index);
1184  }
1185 
1186 
1187  inline bool
1188  const_iterator::operator!=(const const_iterator &other) const
1189  {
1190  return !(*this == other);
1191  }
1192 
1193 
1194  inline bool
1195  const_iterator::operator<(const const_iterator &other) const
1196  {
1197  return (accessor.row() < other.accessor.row() ||
1198  (accessor.row() == other.accessor.row() &&
1199  accessor.index() < other.accessor.index()));
1200  }
1201 
1202  } // namespace MatrixIterators
1203 
1204 
1205 
1206  // Inline the set() and add()
1207  // functions, since they will be
1208  // called frequently, and the
1209  // compiler can optimize away
1210  // some unnecessary loops when
1211  // the sizes are given at
1212  // compile time.
1213  inline void
1214  MatrixBase::set(const size_type i, const size_type j, const PetscScalar value)
1215  {
1216  AssertIsFinite(value);
1217 
1218  set(i, 1, &j, &value, false);
1219  }
1220 
1221 
1222 
1223  inline void
1224  MatrixBase::set(const std::vector<size_type> & indices,
1225  const FullMatrix<PetscScalar> &values,
1226  const bool elide_zero_values)
1227  {
1228  Assert(indices.size() == values.m(),
1229  ExcDimensionMismatch(indices.size(), values.m()));
1230  Assert(values.m() == values.n(), ExcNotQuadratic());
1231 
1232  for (size_type i = 0; i < indices.size(); ++i)
1233  set(indices[i],
1234  indices.size(),
1235  indices.data(),
1236  &values(i, 0),
1237  elide_zero_values);
1238  }
1239 
1240 
1241 
1242  inline void
1243  MatrixBase::set(const std::vector<size_type> & row_indices,
1244  const std::vector<size_type> & col_indices,
1245  const FullMatrix<PetscScalar> &values,
1246  const bool elide_zero_values)
1247  {
1248  Assert(row_indices.size() == values.m(),
1249  ExcDimensionMismatch(row_indices.size(), values.m()));
1250  Assert(col_indices.size() == values.n(),
1251  ExcDimensionMismatch(col_indices.size(), values.n()));
1252 
1253  for (size_type i = 0; i < row_indices.size(); ++i)
1254  set(row_indices[i],
1255  col_indices.size(),
1256  col_indices.data(),
1257  &values(i, 0),
1258  elide_zero_values);
1259  }
1260 
1261 
1262 
1263  inline void
1264  MatrixBase::set(const size_type row,
1265  const std::vector<size_type> & col_indices,
1266  const std::vector<PetscScalar> &values,
1267  const bool elide_zero_values)
1268  {
1269  Assert(col_indices.size() == values.size(),
1270  ExcDimensionMismatch(col_indices.size(), values.size()));
1271 
1272  set(row,
1273  col_indices.size(),
1274  col_indices.data(),
1275  values.data(),
1276  elide_zero_values);
1277  }
1278 
1279 
1280 
1281  inline void
1282  MatrixBase::set(const size_type row,
1283  const size_type n_cols,
1284  const size_type * col_indices,
1285  const PetscScalar *values,
1286  const bool elide_zero_values)
1287  {
1288  prepare_action(VectorOperation::insert);
1289 
1290  const PetscInt petsc_i = row;
1291  PetscInt const *col_index_ptr;
1292 
1293  PetscScalar const *col_value_ptr;
1294  int n_columns;
1295 
1296  // If we don't elide zeros, the pointers are already available...
1297  if (elide_zero_values == false)
1298  {
1299  col_index_ptr = reinterpret_cast<const PetscInt *>(col_indices);
1300  col_value_ptr = values;
1301  n_columns = n_cols;
1302  }
1303  else
1304  {
1305  // Otherwise, extract nonzero values in each row and get the
1306  // respective index.
1307  if (column_indices.size() < n_cols)
1308  {
1309  column_indices.resize(n_cols);
1310  column_values.resize(n_cols);
1311  }
1312 
1313  n_columns = 0;
1314  for (size_type j = 0; j < n_cols; ++j)
1315  {
1316  const PetscScalar value = values[j];
1317  AssertIsFinite(value);
1318  if (value != PetscScalar())
1319  {
1320  column_indices[n_columns] = col_indices[j];
1321  column_values[n_columns] = value;
1322  n_columns++;
1323  }
1324  }
1325  AssertIndexRange(n_columns, n_cols + 1);
1326 
1327  col_index_ptr = column_indices.data();
1328  col_value_ptr = column_values.data();
1329  }
1330 
1331  const PetscErrorCode ierr = MatSetValues(matrix,
1332  1,
1333  &petsc_i,
1334  n_columns,
1335  col_index_ptr,
1336  col_value_ptr,
1337  INSERT_VALUES);
1338  AssertThrow(ierr == 0, ExcPETScError(ierr));
1339  }
1340 
1341 
1342 
1343  inline void
1344  MatrixBase::add(const size_type i, const size_type j, const PetscScalar value)
1345  {
1346  AssertIsFinite(value);
1347 
1348  if (value == PetscScalar())
1349  {
1350  // we have to check after using Insert/Add in any case to be
1351  // consistent with the MPI communication model, but we can save
1352  // some work if the addend is zero. However, these actions are done
1353  // in case we pass on to the other function.
1354  prepare_action(VectorOperation::add);
1355 
1356  return;
1357  }
1358  else
1359  add(i, 1, &j, &value, false);
1360  }
1361 
1362 
1363 
1364  inline void
1365  MatrixBase::add(const std::vector<size_type> & indices,
1366  const FullMatrix<PetscScalar> &values,
1367  const bool elide_zero_values)
1368  {
1369  Assert(indices.size() == values.m(),
1370  ExcDimensionMismatch(indices.size(), values.m()));
1371  Assert(values.m() == values.n(), ExcNotQuadratic());
1372 
1373  for (size_type i = 0; i < indices.size(); ++i)
1374  add(indices[i],
1375  indices.size(),
1376  indices.data(),
1377  &values(i, 0),
1378  elide_zero_values);
1379  }
1380 
1381 
1382 
1383  inline void
1384  MatrixBase::add(const std::vector<size_type> & row_indices,
1385  const std::vector<size_type> & col_indices,
1386  const FullMatrix<PetscScalar> &values,
1387  const bool elide_zero_values)
1388  {
1389  Assert(row_indices.size() == values.m(),
1390  ExcDimensionMismatch(row_indices.size(), values.m()));
1391  Assert(col_indices.size() == values.n(),
1392  ExcDimensionMismatch(col_indices.size(), values.n()));
1393 
1394  for (size_type i = 0; i < row_indices.size(); ++i)
1395  add(row_indices[i],
1396  col_indices.size(),
1397  col_indices.data(),
1398  &values(i, 0),
1399  elide_zero_values);
1400  }
1401 
1402 
1403 
1404  inline void
1405  MatrixBase::add(const size_type row,
1406  const std::vector<size_type> & col_indices,
1407  const std::vector<PetscScalar> &values,
1408  const bool elide_zero_values)
1409  {
1410  Assert(col_indices.size() == values.size(),
1411  ExcDimensionMismatch(col_indices.size(), values.size()));
1412 
1413  add(row,
1414  col_indices.size(),
1415  col_indices.data(),
1416  values.data(),
1417  elide_zero_values);
1418  }
1419 
1420 
1421 
1422  inline void
1423  MatrixBase::add(const size_type row,
1424  const size_type n_cols,
1425  const size_type * col_indices,
1426  const PetscScalar *values,
1427  const bool elide_zero_values,
1428  const bool /*col_indices_are_sorted*/)
1429  {
1430  (void)elide_zero_values;
1431 
1432  prepare_action(VectorOperation::add);
1433 
1434  const PetscInt petsc_i = row;
1435  PetscInt const *col_index_ptr;
1436 
1437  PetscScalar const *col_value_ptr;
1438  int n_columns;
1439 
1440  // If we don't elide zeros, the pointers are already available...
1441  if (elide_zero_values == false)
1442  {
1443  col_index_ptr = reinterpret_cast<const PetscInt *>(col_indices);
1444  col_value_ptr = values;
1445  n_columns = n_cols;
1446  }
1447  else
1448  {
1449  // Otherwise, extract nonzero values in each row and get the
1450  // respective index.
1451  if (column_indices.size() < n_cols)
1452  {
1453  column_indices.resize(n_cols);
1454  column_values.resize(n_cols);
1455  }
1456 
1457  n_columns = 0;
1458  for (size_type j = 0; j < n_cols; ++j)
1459  {
1460  const PetscScalar value = values[j];
1461  AssertIsFinite(value);
1462  if (value != PetscScalar())
1463  {
1464  column_indices[n_columns] = col_indices[j];
1465  column_values[n_columns] = value;
1466  n_columns++;
1467  }
1468  }
1469  AssertIndexRange(n_columns, n_cols + 1);
1470 
1471  col_index_ptr = column_indices.data();
1472  col_value_ptr = column_values.data();
1473  }
1474 
1475  const PetscErrorCode ierr = MatSetValues(
1476  matrix, 1, &petsc_i, n_columns, col_index_ptr, col_value_ptr, ADD_VALUES);
1477  AssertThrow(ierr == 0, ExcPETScError(ierr));
1478  }
1479 
1480 
1481 
1482  inline PetscScalar
1483  MatrixBase::operator()(const size_type i, const size_type j) const
1484  {
1485  return el(i, j);
1486  }
1487 
1488 
1489 
1491  MatrixBase::begin() const
1492  {
1493  Assert(
1494  (in_local_range(0) && in_local_range(m() - 1)),
1495  ExcMessage(
1496  "begin() and end() can only be called on a processor owning the entire matrix. If this is a distributed matrix, use begin(row) and end(row) instead."));
1497 
1498  // find the first non-empty row in order to make sure that the returned
1499  // iterator points to something useful
1500  size_type first_nonempty_row = 0;
1501  while ((first_nonempty_row < m()) && (row_length(first_nonempty_row) == 0))
1502  ++first_nonempty_row;
1503 
1504  return const_iterator(this, first_nonempty_row, 0);
1505  }
1506 
1507 
1509  MatrixBase::end() const
1510  {
1511  Assert(
1512  (in_local_range(0) && in_local_range(m() - 1)),
1513  ExcMessage(
1514  "begin() and end() can only be called on a processor owning the entire matrix. If this is a distributed matrix, use begin(row) and end(row) instead."));
1515 
1516  return const_iterator(this, m(), 0);
1517  }
1518 
1519 
1521  MatrixBase::begin(const size_type r) const
1522  {
1523  Assert(in_local_range(r),
1524  ExcIndexRange(r, local_range().first, local_range().second));
1525 
1526  if (row_length(r) > 0)
1527  return const_iterator(this, r, 0);
1528  else
1529  return end(r);
1530  }
1531 
1532 
1534  MatrixBase::end(const size_type r) const
1535  {
1536  Assert(in_local_range(r),
1537  ExcIndexRange(r, local_range().first, local_range().second));
1538 
1539  // place the iterator on the first entry past this line, or at the
1540  // end of the matrix
1541  //
1542  // in the parallel case, we need to put it on the first entry of
1543  // the first row after the locally owned range. this of course
1544  // doesn't exist, but we can nevertheless create such an
1545  // iterator. we need to check whether 'i' is past the locally
1546  // owned range of rows first, before we ask for the length of the
1547  // row since the latter query leads to an exception in case we ask
1548  // for a row that is not locally owned
1549  for (size_type i = r + 1; i < m(); ++i)
1550  if (i == local_range().second || (row_length(i) > 0))
1551  return const_iterator(this, i, 0);
1552 
1553  // if there is no such line, then take the
1554  // end iterator of the matrix
1555  // we don't allow calling end() directly for distributed matrices so we need
1556  // to copy the code without the assertion.
1557  return {this, m(), 0};
1558  }
1559 
1560 
1561 
1562  inline bool
1564  {
1565  PetscInt begin, end;
1566 
1567  const PetscErrorCode ierr =
1568  MatGetOwnershipRange(static_cast<const Mat &>(matrix), &begin, &end);
1569  AssertThrow(ierr == 0, ExcPETScError(ierr));
1570 
1571  return ((index >= static_cast<size_type>(begin)) &&
1572  (index < static_cast<size_type>(end)));
1573  }
1574 
1575 
1576 
1577  inline void
1579  {
1580  if (last_action == VectorOperation::unknown)
1581  last_action = new_action;
1582 
1583  Assert(last_action == new_action, ExcWrongMode(last_action, new_action));
1584  }
1585 
1586 
1587 
1588  inline void
1590  {
1591  // compress() sets the last action to none, which allows us to check if
1592  // there are pending add/insert operations:
1593  AssertThrow(last_action == VectorOperation::unknown,
1594  ExcMessage("Error: missing compress() call."));
1595  }
1596 
1597 
1598 
1599  inline void
1601  {
1602  prepare_action(VectorOperation::add);
1603  }
1604 
1605 
1606 
1607  inline void
1609  {
1610  prepare_action(VectorOperation::insert);
1611  }
1612 
1613 # endif // DOXYGEN
1614 } // namespace PETScWrappers
1615 
1616 
1618 
1619 
1620 # endif // DEAL_II_WITH_PETSC
1621 
1622 #endif
1623 /*---------------------------- petsc_matrix_base.h --------------------------*/
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:532
bool operator==(const const_iterator &) const
void add(const size_type i, const size_type j, const PetscScalar value)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
const_iterator begin() const
static ::ExceptionBase & ExcInvalidIndexWithinRow(int arg1, int arg2)
std::shared_ptr< const std::vector< PetscScalar > > value_cache
#define AssertIndexRange(index, range)
Definition: exceptions.h:1718
static const char V
size_type row_length(const size_type row) const
#define AssertThrow(cond, exc)
Definition: exceptions.h:1571
Point< 2 > second
Definition: grid_out.cc:4587
SymmetricTensor< 2, dim, Number > C(const Tensor< 2, dim, Number > &F)
std::vector< PetscInt > column_indices
const_iterator(const MatrixBase *matrix, const size_type row, const size_type index)
std::vector< PetscScalar > column_values
static ::ExceptionBase & ExcMessage(std::string arg1)
std::string compress(const std::string &input)
Definition: utilities.cc:392
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:1461
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
Number linfty_norm(const Tensor< 2, dim, Number > &t)
Definition: tensor.h:3025
#define DeclException0(Exception0)
Definition: exceptions.h:464
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
VectorType::value_type * end(VectorType &V)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
const_iterator end() const
PetscScalar operator()(const size_type i, const size_type j) const
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
Point< 2 > first
Definition: grid_out.cc:4586
std::shared_ptr< const std::vector< size_type > > colnum_cache
static ::ExceptionBase & ExcIteratorPastEnd()
static ::ExceptionBase & ExcNotQuadratic()
VectorOperation::values last_action
void prepare_action(const VectorOperation::values new_action)
unsigned int global_dof_index
Definition: types.h:76
bool in_local_range(const size_type index) const
std::pair< size_type, size_type > local_range() const
static ::ExceptionBase & ExcAccessToNonlocalRow(int arg1, int arg2, int arg3)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
VectorType::value_type * begin(VectorType &V)
Number l1_norm(const Tensor< 2, dim, Number > &t)
Definition: tensor.h:2999
void set(const size_type i, const size_type j, const PetscScalar value)
constexpr Number trace(const SymmetricTensor< 2, dim2, Number > &)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition: exceptions.h:555
bool operator!=(const const_iterator &) const
#define AssertIsFinite(number)
Definition: exceptions.h:1744
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)