Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
petsc_matrix_base.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2004 - 2023 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
31
32# include <petscmat.h>
33
34# include <cmath>
35# include <memory>
36# include <vector>
37
39
40// Forward declarations
41# ifndef DOXYGEN
42template <typename Matrix>
43class BlockMatrixBase;
44# endif
45
46
47namespace PETScWrappers
48{
49 // forward declarations
50 class MatrixBase;
51
52 namespace MatrixIterators
53 {
68 {
69 private:
74 {
75 public:
80
86 const size_type row,
87 const size_type index);
88
93 row() const;
94
99 index() const;
100
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:
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
195
201
205 const Accessor &
206 operator*() const;
207
211 const Accessor *
212 operator->() const;
213
218 bool
223 bool
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
325 explicit MatrixBase(const Mat &);
326
332 MatrixBase(const MatrixBase &) = delete;
333
339 MatrixBase &
340 operator=(const MatrixBase &) = delete;
341
345 virtual ~MatrixBase() override;
346
353 void
354 reinit(Mat A);
355
365 MatrixBase &
366 operator=(const value_type d);
367
372 void
373 clear();
374
384 void
385 set(const size_type i, const size_type j, const PetscScalar value);
386
406 void
407 set(const std::vector<size_type> & indices,
408 const FullMatrix<PetscScalar> &full_matrix,
409 const bool elide_zero_values = false);
410
416 void
417 set(const std::vector<size_type> & row_indices,
418 const std::vector<size_type> & col_indices,
419 const FullMatrix<PetscScalar> &full_matrix,
420 const bool elide_zero_values = false);
421
436 void
437 set(const size_type row,
438 const std::vector<size_type> & col_indices,
439 const std::vector<PetscScalar> &values,
440 const bool elide_zero_values = false);
441
456 void
457 set(const size_type row,
458 const size_type n_cols,
459 const size_type * col_indices,
460 const PetscScalar *values,
461 const bool elide_zero_values = false);
462
472 void
473 add(const size_type i, const size_type j, const PetscScalar value);
474
494 void
495 add(const std::vector<size_type> & indices,
496 const FullMatrix<PetscScalar> &full_matrix,
497 const bool elide_zero_values = true);
498
504 void
505 add(const std::vector<size_type> & row_indices,
506 const std::vector<size_type> & col_indices,
507 const FullMatrix<PetscScalar> &full_matrix,
508 const bool elide_zero_values = true);
509
524 void
525 add(const size_type row,
526 const std::vector<size_type> & col_indices,
527 const std::vector<PetscScalar> &values,
528 const bool elide_zero_values = true);
529
544 void
545 add(const size_type row,
546 const size_type n_cols,
547 const size_type * col_indices,
548 const PetscScalar *values,
549 const bool elide_zero_values = true,
550 const bool col_indices_are_sorted = false);
551
568 void
569 clear_row(const size_type row, const PetscScalar new_diag_value = 0);
570
579 void
580 clear_rows(const std::vector<size_type> &rows,
581 const PetscScalar new_diag_value = 0);
582
586 void
587 clear_rows_columns(const std::vector<size_type> &row_and_column_indices,
588 const PetscScalar new_diag_value = 0);
589
601 void
602 compress(const VectorOperation::values operation);
603
615 PetscScalar
616 operator()(const size_type i, const size_type j) const;
617
625 PetscScalar
626 el(const size_type i, const size_type j) const;
627
637 PetscScalar
638 diag_element(const size_type i) const;
639
644 m() const;
645
650 n() const;
651
661 local_size() const;
662
671 std::pair<size_type, size_type>
672 local_range() const;
673
678 bool
679 in_local_range(const size_type index) const;
680
688 local_domain_size() const;
689
696 std::pair<size_type, size_type>
697 local_domain() const;
698
704
710 std::uint64_t
711 n_nonzero_elements() const;
712
717 row_length(const size_type row) const;
718
726 PetscReal
727 l1_norm() const;
728
736 PetscReal
737 linfty_norm() const;
738
743 PetscReal
744 frobenius_norm() const;
745
746
766 PetscScalar
767 matrix_norm_square(const VectorBase &v) const;
768
769
783 PetscScalar
784 matrix_scalar_product(const VectorBase &u, const VectorBase &v) const;
785
790 PetscScalar
791 trace() const;
792
796 MatrixBase &
797 operator*=(const PetscScalar factor);
798
802 MatrixBase &
803 operator/=(const PetscScalar factor);
804
805
810 MatrixBase &
811 add(const PetscScalar factor, const MatrixBase &other);
812
824 void
825 vmult(VectorBase &dst, const VectorBase &src) const;
826
839 void
840 Tvmult(VectorBase &dst, const VectorBase &src) const;
841
853 void
854 vmult_add(VectorBase &dst, const VectorBase &src) const;
855
868 void
869 Tvmult_add(VectorBase &dst, const VectorBase &src) const;
870
883 PetscScalar
884 residual(VectorBase &dst, const VectorBase &x, const VectorBase &b) const;
885
892 begin() const;
893
900 end() const;
901
911 begin(const size_type r) const;
912
922 end(const size_type r) const;
923
931 operator Mat() const;
932
938 Mat &
939 petsc_matrix();
940
944 void
945 transpose();
946
951 PetscBool
952 is_symmetric(const double tolerance = 1.e-12);
953
959 PetscBool
960 is_hermitian(const double tolerance = 1.e-12);
961
968 void
969 write_ascii(const PetscViewerFormat format = PETSC_VIEWER_DEFAULT);
970
978 void
979 print(std::ostream &out, const bool alternative_output = false) const;
980
984 std::size_t
985 memory_consumption() const;
986
991 "You are attempting an operation on two vectors that "
992 "are the same object, but the operation requires that the "
993 "two objects are in fact different.");
994
999 int,
1000 int,
1001 << "You tried to do a "
1002 << (arg1 == 1 ? "'set'" : (arg1 == 2 ? "'add'" : "???"))
1003 << " operation but the matrix is currently in "
1004 << (arg2 == 1 ? "'set'" : (arg2 == 2 ? "'add'" : "???"))
1005 << " mode. You first have to call 'compress()'.");
1006
1007 protected:
1013
1018
1024 void
1026
1032 void
1034
1045 void
1052 void
1054
1070 void
1071 mmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const;
1072
1089 void
1090 Tmmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const;
1091
1092 private:
1107 mutable std::vector<PetscInt> column_indices;
1108
1117 mutable std::vector<PetscScalar> column_values;
1118
1119
1120 // To allow calling protected prepare_add() and prepare_set().
1121 template <class>
1122 friend class ::BlockMatrixBase;
1123 };
1124
1125
1126
1127# ifndef DOXYGEN
1128 // ---------------------- inline and template functions ---------------------
1129
1130
1131 namespace MatrixIterators
1132 {
1134 const size_type row,
1135 const size_type index)
1136 : matrix(const_cast<MatrixBase *>(matrix))
1137 , a_row(row)
1138 , a_index(index)
1139 {
1140 visit_present_row();
1141 }
1142
1143
1144
1147 {
1148 Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1149 return a_row;
1150 }
1151
1152
1155 {
1156 Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1157 return (*colnum_cache)[a_index];
1158 }
1159
1160
1163 {
1164 Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1165 return a_index;
1166 }
1167
1168
1169 inline PetscScalar
1171 {
1172 Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1173 return (*value_cache)[a_index];
1174 }
1175
1176
1177 inline const_iterator::const_iterator(const MatrixBase *matrix,
1178 const size_type row,
1179 const size_type index)
1180 : accessor(matrix, row, index)
1181 {}
1182
1183
1184
1185 inline const_iterator &
1186 const_iterator::operator++()
1187 {
1188 Assert(accessor.a_row < accessor.matrix->m(), ExcIteratorPastEnd());
1189
1190 ++accessor.a_index;
1191
1192 // if at end of line: do one step, then cycle until we find a
1193 // row with a nonzero number of entries
1194 if (accessor.a_index >= accessor.colnum_cache->size())
1195 {
1196 accessor.a_index = 0;
1197 ++accessor.a_row;
1198
1199 while ((accessor.a_row < accessor.matrix->m()) &&
1200 (accessor.a_row < accessor.matrix->local_range().second) &&
1201 (accessor.matrix->row_length(accessor.a_row) == 0))
1202 ++accessor.a_row;
1203
1204 accessor.visit_present_row();
1205 }
1206 return *this;
1207 }
1208
1209
1210 inline const_iterator
1211 const_iterator::operator++(int)
1212 {
1213 const const_iterator old_state = *this;
1214 ++(*this);
1215 return old_state;
1216 }
1217
1218
1219 inline const const_iterator::Accessor &
1220 const_iterator::operator*() const
1221 {
1222 return accessor;
1223 }
1224
1225
1226 inline const const_iterator::Accessor *
1227 const_iterator::operator->() const
1228 {
1229 return &accessor;
1230 }
1231
1232
1233 inline bool
1234 const_iterator::operator==(const const_iterator &other) const
1235 {
1236 return (accessor.a_row == other.accessor.a_row &&
1237 accessor.a_index == other.accessor.a_index);
1238 }
1239
1240
1241 inline bool
1242 const_iterator::operator!=(const const_iterator &other) const
1243 {
1244 return !(*this == other);
1245 }
1246
1247
1248 inline bool
1249 const_iterator::operator<(const const_iterator &other) const
1250 {
1251 return (accessor.row() < other.accessor.row() ||
1252 (accessor.row() == other.accessor.row() &&
1253 accessor.index() < other.accessor.index()));
1254 }
1255
1256 } // namespace MatrixIterators
1257
1258
1259
1260 // Inline the set() and add()
1261 // functions, since they will be
1262 // called frequently, and the
1263 // compiler can optimize away
1264 // some unnecessary loops when
1265 // the sizes are given at
1266 // compile time.
1267 inline void
1268 MatrixBase::set(const size_type i, const size_type j, const PetscScalar value)
1269 {
1270 AssertIsFinite(value);
1271
1272 set(i, 1, &j, &value, false);
1273 }
1274
1275
1276
1277 inline void
1278 MatrixBase::set(const std::vector<size_type> & indices,
1279 const FullMatrix<PetscScalar> &values,
1280 const bool elide_zero_values)
1281 {
1282 Assert(indices.size() == values.m(),
1283 ExcDimensionMismatch(indices.size(), values.m()));
1284 Assert(values.m() == values.n(), ExcNotQuadratic());
1285
1286 for (size_type i = 0; i < indices.size(); ++i)
1287 set(indices[i],
1288 indices.size(),
1289 indices.data(),
1290 &values(i, 0),
1291 elide_zero_values);
1292 }
1293
1294
1295
1296 inline void
1297 MatrixBase::set(const std::vector<size_type> & row_indices,
1298 const std::vector<size_type> & col_indices,
1299 const FullMatrix<PetscScalar> &values,
1300 const bool elide_zero_values)
1301 {
1302 Assert(row_indices.size() == values.m(),
1303 ExcDimensionMismatch(row_indices.size(), values.m()));
1304 Assert(col_indices.size() == values.n(),
1305 ExcDimensionMismatch(col_indices.size(), values.n()));
1306
1307 for (size_type i = 0; i < row_indices.size(); ++i)
1308 set(row_indices[i],
1309 col_indices.size(),
1310 col_indices.data(),
1311 &values(i, 0),
1312 elide_zero_values);
1313 }
1314
1315
1316
1317 inline void
1318 MatrixBase::set(const size_type row,
1319 const std::vector<size_type> & col_indices,
1320 const std::vector<PetscScalar> &values,
1321 const bool elide_zero_values)
1322 {
1323 Assert(col_indices.size() == values.size(),
1324 ExcDimensionMismatch(col_indices.size(), values.size()));
1325
1326 set(row,
1327 col_indices.size(),
1328 col_indices.data(),
1329 values.data(),
1330 elide_zero_values);
1331 }
1332
1333
1334
1335 inline void
1336 MatrixBase::set(const size_type row,
1337 const size_type n_cols,
1338 const size_type * col_indices,
1339 const PetscScalar *values,
1340 const bool elide_zero_values)
1341 {
1342 prepare_action(VectorOperation::insert);
1343
1344 const PetscInt petsc_i = row;
1345 PetscInt const *col_index_ptr;
1346
1347 PetscScalar const *col_value_ptr;
1348 int n_columns;
1349
1350 // If we don't elide zeros, the pointers are already available...
1351 if (elide_zero_values == false)
1352 {
1353 col_index_ptr = reinterpret_cast<const PetscInt *>(col_indices);
1354 col_value_ptr = values;
1355 n_columns = n_cols;
1356 }
1357 else
1358 {
1359 // Otherwise, extract nonzero values in each row and get the
1360 // respective index.
1361 if (column_indices.size() < n_cols)
1362 {
1363 column_indices.resize(n_cols);
1364 column_values.resize(n_cols);
1365 }
1366
1367 n_columns = 0;
1368 for (size_type j = 0; j < n_cols; ++j)
1369 {
1370 const PetscScalar value = values[j];
1371 AssertIsFinite(value);
1372 if (value != PetscScalar())
1373 {
1374 column_indices[n_columns] = col_indices[j];
1375 column_values[n_columns] = value;
1376 n_columns++;
1377 }
1378 }
1379 AssertIndexRange(n_columns, n_cols + 1);
1380
1381 col_index_ptr = column_indices.data();
1382 col_value_ptr = column_values.data();
1383 }
1384
1385 const PetscErrorCode ierr = MatSetValues(matrix,
1386 1,
1387 &petsc_i,
1388 n_columns,
1389 col_index_ptr,
1390 col_value_ptr,
1391 INSERT_VALUES);
1392 AssertThrow(ierr == 0, ExcPETScError(ierr));
1393 }
1394
1395
1396
1397 inline void
1398 MatrixBase::add(const size_type i, const size_type j, const PetscScalar value)
1399 {
1400 AssertIsFinite(value);
1401
1402 if (value == PetscScalar())
1403 {
1404 // we have to check after using Insert/Add in any case to be
1405 // consistent with the MPI communication model, but we can save
1406 // some work if the addend is zero. However, these actions are done
1407 // in case we pass on to the other function.
1408 prepare_action(VectorOperation::add);
1409
1410 return;
1411 }
1412 else
1413 add(i, 1, &j, &value, false);
1414 }
1415
1416
1417
1418 inline void
1419 MatrixBase::add(const std::vector<size_type> & indices,
1420 const FullMatrix<PetscScalar> &values,
1421 const bool elide_zero_values)
1422 {
1423 Assert(indices.size() == values.m(),
1424 ExcDimensionMismatch(indices.size(), values.m()));
1425 Assert(values.m() == values.n(), ExcNotQuadratic());
1426
1427 for (size_type i = 0; i < indices.size(); ++i)
1428 add(indices[i],
1429 indices.size(),
1430 indices.data(),
1431 &values(i, 0),
1432 elide_zero_values);
1433 }
1434
1435
1436
1437 inline void
1438 MatrixBase::add(const std::vector<size_type> & row_indices,
1439 const std::vector<size_type> & col_indices,
1440 const FullMatrix<PetscScalar> &values,
1441 const bool elide_zero_values)
1442 {
1443 Assert(row_indices.size() == values.m(),
1444 ExcDimensionMismatch(row_indices.size(), values.m()));
1445 Assert(col_indices.size() == values.n(),
1446 ExcDimensionMismatch(col_indices.size(), values.n()));
1447
1448 for (size_type i = 0; i < row_indices.size(); ++i)
1449 add(row_indices[i],
1450 col_indices.size(),
1451 col_indices.data(),
1452 &values(i, 0),
1453 elide_zero_values);
1454 }
1455
1456
1457
1458 inline void
1459 MatrixBase::add(const size_type row,
1460 const std::vector<size_type> & col_indices,
1461 const std::vector<PetscScalar> &values,
1462 const bool elide_zero_values)
1463 {
1464 Assert(col_indices.size() == values.size(),
1465 ExcDimensionMismatch(col_indices.size(), values.size()));
1466
1467 add(row,
1468 col_indices.size(),
1469 col_indices.data(),
1470 values.data(),
1471 elide_zero_values);
1472 }
1473
1474
1475
1476 inline void
1477 MatrixBase::add(const size_type row,
1478 const size_type n_cols,
1479 const size_type * col_indices,
1480 const PetscScalar *values,
1481 const bool elide_zero_values,
1482 const bool /*col_indices_are_sorted*/)
1483 {
1484 (void)elide_zero_values;
1485
1486 prepare_action(VectorOperation::add);
1487
1488 const PetscInt petsc_i = row;
1489 PetscInt const *col_index_ptr;
1490
1491 PetscScalar const *col_value_ptr;
1492 int n_columns;
1493
1494 // If we don't elide zeros, the pointers are already available...
1495 if (elide_zero_values == false)
1496 {
1497 col_index_ptr = reinterpret_cast<const PetscInt *>(col_indices);
1498 col_value_ptr = values;
1499 n_columns = n_cols;
1500 }
1501 else
1502 {
1503 // Otherwise, extract nonzero values in each row and get the
1504 // respective index.
1505 if (column_indices.size() < n_cols)
1506 {
1507 column_indices.resize(n_cols);
1508 column_values.resize(n_cols);
1509 }
1510
1511 n_columns = 0;
1512 for (size_type j = 0; j < n_cols; ++j)
1513 {
1514 const PetscScalar value = values[j];
1515 AssertIsFinite(value);
1516 if (value != PetscScalar())
1517 {
1518 column_indices[n_columns] = col_indices[j];
1519 column_values[n_columns] = value;
1520 n_columns++;
1521 }
1522 }
1523 AssertIndexRange(n_columns, n_cols + 1);
1524
1525 col_index_ptr = column_indices.data();
1526 col_value_ptr = column_values.data();
1527 }
1528
1529 const PetscErrorCode ierr = MatSetValues(
1530 matrix, 1, &petsc_i, n_columns, col_index_ptr, col_value_ptr, ADD_VALUES);
1531 AssertThrow(ierr == 0, ExcPETScError(ierr));
1532 }
1533
1534
1535
1536 inline PetscScalar
1537 MatrixBase::operator()(const size_type i, const size_type j) const
1538 {
1539 return el(i, j);
1540 }
1541
1542
1543
1544 inline MatrixBase::const_iterator
1545 MatrixBase::begin() const
1546 {
1547 Assert(
1548 (in_local_range(0) && in_local_range(m() - 1)),
1549 ExcMessage(
1550 "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."));
1551
1552 // find the first non-empty row in order to make sure that the returned
1553 // iterator points to something useful
1554 size_type first_nonempty_row = 0;
1555 while ((first_nonempty_row < m()) && (row_length(first_nonempty_row) == 0))
1556 ++first_nonempty_row;
1557
1558 return const_iterator(this, first_nonempty_row, 0);
1559 }
1560
1561
1562 inline MatrixBase::const_iterator
1563 MatrixBase::end() const
1564 {
1565 Assert(
1566 (in_local_range(0) && in_local_range(m() - 1)),
1567 ExcMessage(
1568 "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."));
1569
1570 return const_iterator(this, m(), 0);
1571 }
1572
1573
1574 inline MatrixBase::const_iterator
1575 MatrixBase::begin(const size_type r) const
1576 {
1577 Assert(in_local_range(r),
1578 ExcIndexRange(r, local_range().first, local_range().second));
1579
1580 if (row_length(r) > 0)
1581 return const_iterator(this, r, 0);
1582 else
1583 return end(r);
1584 }
1585
1586
1587 inline MatrixBase::const_iterator
1588 MatrixBase::end(const size_type r) const
1589 {
1590 Assert(in_local_range(r),
1591 ExcIndexRange(r, local_range().first, local_range().second));
1592
1593 // place the iterator on the first entry past this line, or at the
1594 // end of the matrix
1595 //
1596 // in the parallel case, we need to put it on the first entry of
1597 // the first row after the locally owned range. this of course
1598 // doesn't exist, but we can nevertheless create such an
1599 // iterator. we need to check whether 'i' is past the locally
1600 // owned range of rows first, before we ask for the length of the
1601 // row since the latter query leads to an exception in case we ask
1602 // for a row that is not locally owned
1603 for (size_type i = r + 1; i < m(); ++i)
1604 if (i == local_range().second || (row_length(i) > 0))
1605 return const_iterator(this, i, 0);
1606
1607 // if there is no such line, then take the
1608 // end iterator of the matrix
1609 // we don't allow calling end() directly for distributed matrices so we need
1610 // to copy the code without the assertion.
1611 return {this, m(), 0};
1612 }
1613
1614
1615
1616 inline bool
1617 MatrixBase::in_local_range(const size_type index) const
1618 {
1619 PetscInt begin, end;
1620
1621 const PetscErrorCode ierr =
1622 MatGetOwnershipRange(static_cast<const Mat &>(matrix), &begin, &end);
1623 AssertThrow(ierr == 0, ExcPETScError(ierr));
1624
1625 return ((index >= static_cast<size_type>(begin)) &&
1626 (index < static_cast<size_type>(end)));
1627 }
1628
1629
1630
1631 inline void
1632 MatrixBase::prepare_action(const VectorOperation::values new_action)
1633 {
1634 if (last_action == VectorOperation::unknown)
1635 last_action = new_action;
1636
1637 Assert(last_action == new_action, ExcWrongMode(last_action, new_action));
1638 }
1639
1640
1641
1642 inline void
1643 MatrixBase::assert_is_compressed()
1644 {
1645 // compress() sets the last action to none, which allows us to check if
1646 // there are pending add/insert operations:
1648 ExcMessage("Error: missing compress() call."));
1649 }
1650
1651
1652
1653 inline void
1654 MatrixBase::prepare_add()
1655 {
1656 prepare_action(VectorOperation::add);
1657 }
1658
1659
1660
1661 inline void
1662 MatrixBase::prepare_set()
1663 {
1664 prepare_action(VectorOperation::insert);
1665 }
1666
1667 inline MPI_Comm
1668 MatrixBase::get_mpi_communicator() const
1669 {
1670 return PetscObjectComm(reinterpret_cast<PetscObject>(matrix));
1671 }
1672
1673# endif // DOXYGEN
1674} // namespace PETScWrappers
1675
1676
1678
1679
1680#endif // DEAL_II_WITH_PETSC
1681
1682#endif
void add(const size_type i, const size_type j, const PetscScalar value)
void add(const std::vector< size_type > &row_indices, const std::vector< size_type > &col_indices, const FullMatrix< PetscScalar > &full_matrix, const bool elide_zero_values=true)
const_iterator end(const size_type r) const
size_type row_length(const size_type row) const
std::size_t memory_consumption() const
VectorOperation::values last_action
void vmult(VectorBase &dst, const VectorBase &src) const
void clear_rows(const std::vector< size_type > &rows, const PetscScalar new_diag_value=0)
std::vector< PetscScalar > column_values
MPI_Comm get_mpi_communicator() const
PetscScalar diag_element(const size_type i) const
MatrixBase(const MatrixBase &)=delete
size_type local_domain_size() const
const_iterator begin() const
MatrixBase & operator/=(const PetscScalar factor)
void mmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const
PetscBool is_symmetric(const double tolerance=1.e-12)
std::pair< size_type, size_type > local_domain() const
const_iterator end() const
void Tvmult_add(VectorBase &dst, const VectorBase &src) const
void print(std::ostream &out, const bool alternative_output=false) const
const_iterator begin(const size_type r) const
void add(const size_type row, const std::vector< size_type > &col_indices, const std::vector< PetscScalar > &values, const bool elide_zero_values=true)
PetscScalar el(const size_type i, const size_type j) const
PetscScalar operator()(const size_type i, const size_type j) const
MatrixBase & operator=(const MatrixBase &)=delete
void set(const std::vector< size_type > &row_indices, const std::vector< size_type > &col_indices, const FullMatrix< PetscScalar > &full_matrix, const bool elide_zero_values=false)
void prepare_action(const VectorOperation::values new_action)
void add(const std::vector< size_type > &indices, const FullMatrix< PetscScalar > &full_matrix, const bool elide_zero_values=true)
bool in_local_range(const size_type index) const
PetscBool is_hermitian(const double tolerance=1.e-12)
PetscScalar matrix_scalar_product(const VectorBase &u, const VectorBase &v) const
std::vector< PetscInt > column_indices
void Tvmult(VectorBase &dst, const VectorBase &src) const
void clear_rows_columns(const std::vector< size_type > &row_and_column_indices, const PetscScalar new_diag_value=0)
MatrixBase & operator*=(const PetscScalar factor)
PetscScalar residual(VectorBase &dst, const VectorBase &x, const VectorBase &b) const
void Tmmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const
void set(const size_type row, const std::vector< size_type > &col_indices, const std::vector< PetscScalar > &values, const bool elide_zero_values=false)
void add(const size_type row, const size_type n_cols, const size_type *col_indices, const PetscScalar *values, const bool elide_zero_values=true, const bool col_indices_are_sorted=false)
void write_ascii(const PetscViewerFormat format=PETSC_VIEWER_DEFAULT)
void vmult_add(VectorBase &dst, const VectorBase &src) const
std::pair< size_type, size_type > local_range() const
void compress(const VectorOperation::values operation)
PetscScalar matrix_norm_square(const VectorBase &v) const
void set(const size_type row, const size_type n_cols, const size_type *col_indices, const PetscScalar *values, const bool elide_zero_values=false)
void set(const std::vector< size_type > &indices, const FullMatrix< PetscScalar > &full_matrix, const bool elide_zero_values=false)
void set(const size_type i, const size_type j, const PetscScalar value)
std::uint64_t n_nonzero_elements() const
void clear_row(const size_type row, const PetscScalar new_diag_value=0)
std::shared_ptr< const std::vector< PetscScalar > > value_cache
std::shared_ptr< const std::vector< size_type > > colnum_cache
Accessor(const MatrixBase *matrix, const size_type row, const size_type index)
bool operator!=(const const_iterator &) const
const_iterator(const MatrixBase *matrix, const size_type row, const size_type index)
bool operator==(const const_iterator &) const
bool operator<(const const_iterator &) const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
Point< 2 > second
Definition grid_out.cc:4616
Point< 2 > first
Definition grid_out.cc:4615
__global__ void set(Number *val, const Number s, const size_type N)
#define DeclException0(Exception0)
Definition exceptions.h:465
static ::ExceptionBase & ExcAccessToNonlocalRow(int arg1, int arg2, int arg3)
#define Assert(cond, exc)
static ::ExceptionBase & ExcIteratorPastEnd()
static ::ExceptionBase & ExcInvalidIndexWithinRow(int arg1, int arg2)
#define AssertIsFinite(number)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition exceptions.h:533
#define AssertIndexRange(index, range)
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:488
static ::ExceptionBase & ExcSourceEqualsDestination()
static ::ExceptionBase & ExcWrongMode(int arg1, int arg2)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition exceptions.h:556
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcNotQuadratic()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
@ matrix
Contents is actually a matrix.
types::global_dof_index size_type
VectorType::value_type * end(VectorType &V)
VectorType::value_type * begin(VectorType &V)
unsigned int global_dof_index
Definition types.h:82