Reference documentation for deal.II version GIT relicensing-43-g3fc2fae7c3 2024-03-03 21:10:02+00:00
\(\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// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2004 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_petsc_matrix_base_h
16#define dealii_petsc_matrix_base_h
17
18
19#include <deal.II/base/config.h>
20
21#ifdef DEAL_II_WITH_PETSC
22
24
30
31# include <petscmat.h>
32
33# include <cmath>
34# include <memory>
35# include <vector>
36
38
39// Forward declarations
40# ifndef DOXYGEN
41template <typename Matrix>
42class BlockMatrixBase;
43# endif
44
45
46namespace PETScWrappers
47{
48 // forward declarations
49 class MatrixBase;
50
51 namespace MatrixIterators
52 {
67 {
68 private:
73 {
74 public:
79
85 const size_type row,
86 const size_type index);
87
92 row() const;
93
98 index() const;
99
104 column() const;
105
109 PetscScalar
110 value() const;
111
120 int,
121 int,
122 int,
123 << "You tried to access row " << arg1
124 << " of a distributed matrix, but only rows " << arg2
125 << " through " << arg3
126 << " are stored locally and can be accessed.");
127
128 private:
133
138
143
156 std::shared_ptr<const std::vector<size_type>> colnum_cache;
157
161 std::shared_ptr<const std::vector<PetscScalar>> value_cache;
162
168 void
170
171 // Make enclosing class a friend.
172 friend class const_iterator;
173 };
174
175 public:
180
186 const size_type row,
187 const size_type index);
188
194
200
204 const Accessor &
205 operator*() const;
206
210 const Accessor *
211 operator->() const;
212
217 bool
222 bool
224
230 bool
231 operator<(const const_iterator &) const;
232
237 int,
238 int,
239 << "Attempt to access element " << arg2 << " of row "
240 << arg1 << " which doesn't have that many elements.");
241
242 private:
247 };
248
249 } // namespace MatrixIterators
250
251
283 class MatrixBase : public Subscriptor
284 {
285 public:
290
295
299 using value_type = PetscScalar;
300
304 MatrixBase();
305
324 explicit MatrixBase(const Mat &);
325
331 MatrixBase(const MatrixBase &) = delete;
332
338 MatrixBase &
339 operator=(const MatrixBase &) = delete;
340
344 virtual ~MatrixBase() override;
345
352 void
353 reinit(Mat A);
354
364 MatrixBase &
365 operator=(const value_type d);
366
371 void
372 clear();
373
383 void
384 set(const size_type i, const size_type j, const PetscScalar value);
385
405 void
406 set(const std::vector<size_type> &indices,
407 const FullMatrix<PetscScalar> &full_matrix,
408 const bool elide_zero_values = false);
409
415 void
416 set(const std::vector<size_type> &row_indices,
417 const std::vector<size_type> &col_indices,
418 const FullMatrix<PetscScalar> &full_matrix,
419 const bool elide_zero_values = false);
420
435 void
436 set(const size_type row,
437 const std::vector<size_type> &col_indices,
438 const std::vector<PetscScalar> &values,
439 const bool elide_zero_values = false);
440
455 void
456 set(const size_type row,
457 const size_type n_cols,
458 const size_type *col_indices,
459 const PetscScalar *values,
460 const bool elide_zero_values = false);
461
471 void
472 add(const size_type i, const size_type j, const PetscScalar value);
473
493 void
494 add(const std::vector<size_type> &indices,
495 const FullMatrix<PetscScalar> &full_matrix,
496 const bool elide_zero_values = true);
497
503 void
504 add(const std::vector<size_type> &row_indices,
505 const std::vector<size_type> &col_indices,
506 const FullMatrix<PetscScalar> &full_matrix,
507 const bool elide_zero_values = true);
508
523 void
524 add(const size_type row,
525 const std::vector<size_type> &col_indices,
526 const std::vector<PetscScalar> &values,
527 const bool elide_zero_values = true);
528
543 void
544 add(const size_type row,
545 const size_type n_cols,
546 const size_type *col_indices,
547 const PetscScalar *values,
548 const bool elide_zero_values = true,
549 const bool col_indices_are_sorted = false);
550
567 void
568 clear_row(const size_type row, const PetscScalar new_diag_value = 0);
569
578 void
580 const PetscScalar new_diag_value = 0);
581
585 void
586 clear_rows_columns(const std::vector<size_type> &row_and_column_indices,
587 const PetscScalar new_diag_value = 0);
588
600 void
601 compress(const VectorOperation::values operation);
602
614 PetscScalar
615 operator()(const size_type i, const size_type j) const;
616
624 PetscScalar
625 el(const size_type i, const size_type j) const;
626
636 PetscScalar
637 diag_element(const size_type i) const;
638
643 m() const;
644
649 n() const;
650
660 local_size() const;
661
670 std::pair<size_type, size_type>
671 local_range() const;
672
677 bool
678 in_local_range(const size_type index) const;
679
687 local_domain_size() const;
688
695 std::pair<size_type, size_type>
696 local_domain() const;
697
703
709 std::uint64_t
710 n_nonzero_elements() const;
711
716 row_length(const size_type row) const;
717
725 PetscReal
726 l1_norm() const;
727
735 PetscReal
736 linfty_norm() const;
737
742 PetscReal
743 frobenius_norm() const;
744
745
765 PetscScalar
766 matrix_norm_square(const VectorBase &v) const;
767
768
782 PetscScalar
783 matrix_scalar_product(const VectorBase &u, const VectorBase &v) const;
784
789 PetscScalar
790 trace() const;
791
795 MatrixBase &
796 operator*=(const PetscScalar factor);
797
801 MatrixBase &
802 operator/=(const PetscScalar factor);
803
804
809 MatrixBase &
810 add(const PetscScalar factor, const MatrixBase &other);
811
823 void
824 vmult(VectorBase &dst, const VectorBase &src) const;
825
838 void
839 Tvmult(VectorBase &dst, const VectorBase &src) const;
840
852 void
853 vmult_add(VectorBase &dst, const VectorBase &src) const;
854
867 void
868 Tvmult_add(VectorBase &dst, const VectorBase &src) const;
869
882 PetscScalar
883 residual(VectorBase &dst, const VectorBase &x, const VectorBase &b) const;
884
891 begin() const;
892
899 end() const;
900
910 begin(const size_type r) const;
911
921 end(const size_type r) const;
922
930 operator Mat() const;
931
937 Mat &
938 petsc_matrix();
939
943 void
944 transpose();
945
950 PetscBool
951 is_symmetric(const double tolerance = 1.e-12);
952
958 PetscBool
959 is_hermitian(const double tolerance = 1.e-12);
960
967 void
968 write_ascii(const PetscViewerFormat format = PETSC_VIEWER_DEFAULT);
969
977 void
978 print(std::ostream &out, const bool alternative_output = false) const;
979
983 std::size_t
984 memory_consumption() const;
985
990 "You are attempting an operation on two vectors that "
991 "are the same object, but the operation requires that the "
992 "two objects are in fact different.");
993
998 int,
999 int,
1000 << "You tried to do a "
1001 << (arg1 == 1 ? "'set'" : (arg1 == 2 ? "'add'" : "???"))
1002 << " operation but the matrix is currently in "
1003 << (arg2 == 1 ? "'set'" : (arg2 == 2 ? "'add'" : "???"))
1004 << " mode. You first have to call 'compress()'.");
1005
1006 protected:
1012
1017
1023 void
1025
1031 void
1033
1044 void
1051 void
1053
1069 void
1070 mmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const;
1071
1088 void
1089 Tmmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const;
1090
1091 private:
1106 mutable std::vector<PetscInt> column_indices;
1107
1116 mutable std::vector<PetscScalar> column_values;
1117
1118
1119 // To allow calling protected prepare_add() and prepare_set().
1120 template <class>
1121 friend class ::BlockMatrixBase;
1122 };
1123
1124
1125
1126# ifndef DOXYGEN
1127 // ---------------------- inline and template functions ---------------------
1128
1129
1130 namespace MatrixIterators
1131 {
1133 const size_type row,
1134 const size_type index)
1135 : matrix(const_cast<MatrixBase *>(matrix))
1136 , a_row(row)
1137 , a_index(index)
1138 {
1139 visit_present_row();
1140 }
1141
1142
1143
1146 {
1147 Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1148 return a_row;
1149 }
1150
1151
1154 {
1155 Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1156 return (*colnum_cache)[a_index];
1157 }
1158
1159
1162 {
1163 Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1164 return a_index;
1165 }
1166
1167
1168 inline PetscScalar
1170 {
1171 Assert(a_row < matrix->m(), ExcBeyondEndOfMatrix());
1172 return (*value_cache)[a_index];
1173 }
1174
1175
1176 inline const_iterator::const_iterator(const MatrixBase *matrix,
1177 const size_type row,
1178 const size_type index)
1179 : accessor(matrix, row, index)
1180 {}
1181
1182
1183
1184 inline const_iterator &
1185 const_iterator::operator++()
1186 {
1187 Assert(accessor.a_row < accessor.matrix->m(), ExcIteratorPastEnd());
1188
1189 ++accessor.a_index;
1190
1191 // if at end of line: do one step, then cycle until we find a
1192 // row with a nonzero number of entries
1193 if (accessor.a_index >= accessor.colnum_cache->size())
1194 {
1195 accessor.a_index = 0;
1196 ++accessor.a_row;
1197
1198 while ((accessor.a_row < accessor.matrix->m()) &&
1199 (accessor.a_row < accessor.matrix->local_range().second) &&
1200 (accessor.matrix->row_length(accessor.a_row) == 0))
1201 ++accessor.a_row;
1202
1203 accessor.visit_present_row();
1204 }
1205 return *this;
1206 }
1207
1208
1209 inline const_iterator
1210 const_iterator::operator++(int)
1211 {
1212 const const_iterator old_state = *this;
1213 ++(*this);
1214 return old_state;
1215 }
1216
1217
1218 inline const const_iterator::Accessor &
1219 const_iterator::operator*() const
1220 {
1221 return accessor;
1222 }
1223
1224
1225 inline const const_iterator::Accessor *
1226 const_iterator::operator->() const
1227 {
1228 return &accessor;
1229 }
1230
1231
1232 inline bool
1233 const_iterator::operator==(const const_iterator &other) const
1234 {
1235 return (accessor.a_row == other.accessor.a_row &&
1236 accessor.a_index == other.accessor.a_index);
1237 }
1238
1239
1240 inline bool
1241 const_iterator::operator!=(const const_iterator &other) const
1242 {
1243 return !(*this == other);
1244 }
1245
1246
1247 inline bool
1248 const_iterator::operator<(const const_iterator &other) const
1249 {
1250 return (accessor.row() < other.accessor.row() ||
1251 (accessor.row() == other.accessor.row() &&
1252 accessor.index() < other.accessor.index()));
1253 }
1254
1255 } // namespace MatrixIterators
1256
1257
1258
1259 // Inline the set() and add()
1260 // functions, since they will be
1261 // called frequently, and the
1262 // compiler can optimize away
1263 // some unnecessary loops when
1264 // the sizes are given at
1265 // compile time.
1266 inline void
1267 MatrixBase::set(const size_type i, const size_type j, const PetscScalar value)
1268 {
1269 AssertIsFinite(value);
1270
1271 set(i, 1, &j, &value, false);
1272 }
1273
1274
1275
1276 inline void
1277 MatrixBase::set(const std::vector<size_type> &indices,
1278 const FullMatrix<PetscScalar> &values,
1279 const bool elide_zero_values)
1280 {
1281 Assert(indices.size() == values.m(),
1282 ExcDimensionMismatch(indices.size(), values.m()));
1283 Assert(values.m() == values.n(), ExcNotQuadratic());
1284
1285 for (size_type i = 0; i < indices.size(); ++i)
1286 set(indices[i],
1287 indices.size(),
1288 indices.data(),
1289 &values(i, 0),
1290 elide_zero_values);
1291 }
1292
1293
1294
1295 inline void
1296 MatrixBase::set(const std::vector<size_type> &row_indices,
1297 const std::vector<size_type> &col_indices,
1298 const FullMatrix<PetscScalar> &values,
1299 const bool elide_zero_values)
1300 {
1301 Assert(row_indices.size() == values.m(),
1302 ExcDimensionMismatch(row_indices.size(), values.m()));
1303 Assert(col_indices.size() == values.n(),
1304 ExcDimensionMismatch(col_indices.size(), values.n()));
1305
1306 for (size_type i = 0; i < row_indices.size(); ++i)
1307 set(row_indices[i],
1308 col_indices.size(),
1309 col_indices.data(),
1310 &values(i, 0),
1311 elide_zero_values);
1312 }
1313
1314
1315
1316 inline void
1317 MatrixBase::set(const size_type row,
1318 const std::vector<size_type> &col_indices,
1319 const std::vector<PetscScalar> &values,
1320 const bool elide_zero_values)
1321 {
1322 Assert(col_indices.size() == values.size(),
1323 ExcDimensionMismatch(col_indices.size(), values.size()));
1324
1325 set(row,
1326 col_indices.size(),
1327 col_indices.data(),
1328 values.data(),
1329 elide_zero_values);
1330 }
1331
1332
1333
1334 inline void
1335 MatrixBase::set(const size_type row,
1336 const size_type n_cols,
1337 const size_type *col_indices,
1338 const PetscScalar *values,
1339 const bool elide_zero_values)
1340 {
1341 prepare_action(VectorOperation::insert);
1342
1343 const PetscInt petsc_i = row;
1344 const PetscInt *col_index_ptr;
1345
1346 const PetscScalar *col_value_ptr;
1347 int n_columns;
1348
1349 // If we don't elide zeros, the pointers are already available...
1350 if (elide_zero_values == false)
1351 {
1352 col_index_ptr = reinterpret_cast<const PetscInt *>(col_indices);
1353 col_value_ptr = values;
1354 n_columns = n_cols;
1355 }
1356 else
1357 {
1358 // Otherwise, extract nonzero values in each row and get the
1359 // respective index.
1360 if (column_indices.size() < n_cols)
1361 {
1362 column_indices.resize(n_cols);
1363 column_values.resize(n_cols);
1364 }
1365
1366 n_columns = 0;
1367 for (size_type j = 0; j < n_cols; ++j)
1368 {
1369 const PetscScalar value = values[j];
1370 AssertIsFinite(value);
1371 if (value != PetscScalar())
1372 {
1373 column_indices[n_columns] = col_indices[j];
1374 column_values[n_columns] = value;
1375 ++n_columns;
1376 }
1377 }
1378 AssertIndexRange(n_columns, n_cols + 1);
1379
1380 col_index_ptr = column_indices.data();
1381 col_value_ptr = column_values.data();
1382 }
1383
1384 const PetscErrorCode ierr = MatSetValues(matrix,
1385 1,
1386 &petsc_i,
1387 n_columns,
1388 col_index_ptr,
1389 col_value_ptr,
1390 INSERT_VALUES);
1391 AssertThrow(ierr == 0, ExcPETScError(ierr));
1392 }
1393
1394
1395
1396 inline void
1397 MatrixBase::add(const size_type i, const size_type j, const PetscScalar value)
1398 {
1399 AssertIsFinite(value);
1400
1401 if (value == PetscScalar())
1402 {
1403 // we have to check after using Insert/Add in any case to be
1404 // consistent with the MPI communication model, but we can save
1405 // some work if the addend is zero. However, these actions are done
1406 // in case we pass on to the other function.
1407 prepare_action(VectorOperation::add);
1408
1409 return;
1410 }
1411 else
1412 add(i, 1, &j, &value, false);
1413 }
1414
1415
1416
1417 inline void
1418 MatrixBase::add(const std::vector<size_type> &indices,
1419 const FullMatrix<PetscScalar> &values,
1420 const bool elide_zero_values)
1421 {
1422 Assert(indices.size() == values.m(),
1423 ExcDimensionMismatch(indices.size(), values.m()));
1424 Assert(values.m() == values.n(), ExcNotQuadratic());
1425
1426 for (size_type i = 0; i < indices.size(); ++i)
1427 add(indices[i],
1428 indices.size(),
1429 indices.data(),
1430 &values(i, 0),
1431 elide_zero_values);
1432 }
1433
1434
1435
1436 inline void
1437 MatrixBase::add(const std::vector<size_type> &row_indices,
1438 const std::vector<size_type> &col_indices,
1439 const FullMatrix<PetscScalar> &values,
1440 const bool elide_zero_values)
1441 {
1442 Assert(row_indices.size() == values.m(),
1443 ExcDimensionMismatch(row_indices.size(), values.m()));
1444 Assert(col_indices.size() == values.n(),
1445 ExcDimensionMismatch(col_indices.size(), values.n()));
1446
1447 for (size_type i = 0; i < row_indices.size(); ++i)
1448 add(row_indices[i],
1449 col_indices.size(),
1450 col_indices.data(),
1451 &values(i, 0),
1452 elide_zero_values);
1453 }
1454
1455
1456
1457 inline void
1458 MatrixBase::add(const size_type row,
1459 const std::vector<size_type> &col_indices,
1460 const std::vector<PetscScalar> &values,
1461 const bool elide_zero_values)
1462 {
1463 Assert(col_indices.size() == values.size(),
1464 ExcDimensionMismatch(col_indices.size(), values.size()));
1465
1466 add(row,
1467 col_indices.size(),
1468 col_indices.data(),
1469 values.data(),
1470 elide_zero_values);
1471 }
1472
1473
1474
1475 inline void
1476 MatrixBase::add(const size_type row,
1477 const size_type n_cols,
1478 const size_type *col_indices,
1479 const PetscScalar *values,
1480 const bool elide_zero_values,
1481 const bool /*col_indices_are_sorted*/)
1482 {
1483 (void)elide_zero_values;
1484
1485 prepare_action(VectorOperation::add);
1486
1487 const PetscInt petsc_i = row;
1488 const PetscInt *col_index_ptr;
1489
1490 const PetscScalar *col_value_ptr;
1491 int n_columns;
1492
1493 // If we don't elide zeros, the pointers are already available...
1494 if (elide_zero_values == false)
1495 {
1496 col_index_ptr = reinterpret_cast<const PetscInt *>(col_indices);
1497 col_value_ptr = values;
1498 n_columns = n_cols;
1499 }
1500 else
1501 {
1502 // Otherwise, extract nonzero values in each row and get the
1503 // respective index.
1504 if (column_indices.size() < n_cols)
1505 {
1506 column_indices.resize(n_cols);
1507 column_values.resize(n_cols);
1508 }
1509
1510 n_columns = 0;
1511 for (size_type j = 0; j < n_cols; ++j)
1512 {
1513 const PetscScalar value = values[j];
1514 AssertIsFinite(value);
1515 if (value != PetscScalar())
1516 {
1517 column_indices[n_columns] = col_indices[j];
1518 column_values[n_columns] = value;
1519 ++n_columns;
1520 }
1521 }
1522 AssertIndexRange(n_columns, n_cols + 1);
1523
1524 col_index_ptr = column_indices.data();
1525 col_value_ptr = column_values.data();
1526 }
1527
1528 const PetscErrorCode ierr = MatSetValues(
1529 matrix, 1, &petsc_i, n_columns, col_index_ptr, col_value_ptr, ADD_VALUES);
1530 AssertThrow(ierr == 0, ExcPETScError(ierr));
1531 }
1532
1533
1534
1535 inline PetscScalar
1536 MatrixBase::operator()(const size_type i, const size_type j) const
1537 {
1538 return el(i, j);
1539 }
1540
1541
1542
1543 inline MatrixBase::const_iterator
1544 MatrixBase::begin() const
1545 {
1546 Assert(
1547 (in_local_range(0) && in_local_range(m() - 1)),
1548 ExcMessage(
1549 "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."));
1550
1551 // find the first non-empty row in order to make sure that the returned
1552 // iterator points to something useful
1553 size_type first_nonempty_row = 0;
1554 while ((first_nonempty_row < m()) && (row_length(first_nonempty_row) == 0))
1555 ++first_nonempty_row;
1556
1557 return const_iterator(this, first_nonempty_row, 0);
1558 }
1559
1560
1561 inline MatrixBase::const_iterator
1562 MatrixBase::end() const
1563 {
1564 Assert(
1565 (in_local_range(0) && in_local_range(m() - 1)),
1566 ExcMessage(
1567 "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."));
1568
1569 return const_iterator(this, m(), 0);
1570 }
1571
1572
1573 inline MatrixBase::const_iterator
1574 MatrixBase::begin(const size_type r) const
1575 {
1576 Assert(in_local_range(r),
1577 ExcIndexRange(r, local_range().first, local_range().second));
1578
1579 if (row_length(r) > 0)
1580 return const_iterator(this, r, 0);
1581 else
1582 return end(r);
1583 }
1584
1585
1586 inline MatrixBase::const_iterator
1587 MatrixBase::end(const size_type r) const
1588 {
1589 Assert(in_local_range(r),
1590 ExcIndexRange(r, local_range().first, local_range().second));
1591
1592 // place the iterator on the first entry past this line, or at the
1593 // end of the matrix
1594 //
1595 // in the parallel case, we need to put it on the first entry of
1596 // the first row after the locally owned range. this of course
1597 // doesn't exist, but we can nevertheless create such an
1598 // iterator. we need to check whether 'i' is past the locally
1599 // owned range of rows first, before we ask for the length of the
1600 // row since the latter query leads to an exception in case we ask
1601 // for a row that is not locally owned
1602 for (size_type i = r + 1; i < m(); ++i)
1603 if (i == local_range().second || (row_length(i) > 0))
1604 return const_iterator(this, i, 0);
1605
1606 // if there is no such line, then take the
1607 // end iterator of the matrix
1608 // we don't allow calling end() directly for distributed matrices so we need
1609 // to copy the code without the assertion.
1610 return {this, m(), 0};
1611 }
1612
1613
1614
1615 inline bool
1616 MatrixBase::in_local_range(const size_type index) const
1617 {
1618 PetscInt begin, end;
1619
1620 const PetscErrorCode ierr =
1621 MatGetOwnershipRange(static_cast<const Mat &>(matrix), &begin, &end);
1622 AssertThrow(ierr == 0, ExcPETScError(ierr));
1623
1624 return ((index >= static_cast<size_type>(begin)) &&
1625 (index < static_cast<size_type>(end)));
1626 }
1627
1628
1629
1630 inline void
1631 MatrixBase::prepare_action(const VectorOperation::values new_action)
1632 {
1633 if (last_action == VectorOperation::unknown)
1634 last_action = new_action;
1635
1636 Assert(last_action == new_action, ExcWrongMode(last_action, new_action));
1637 }
1638
1639
1640
1641 inline void
1642 MatrixBase::assert_is_compressed()
1643 {
1644 // compress() sets the last action to none, which allows us to check if
1645 // there are pending add/insert operations:
1647 ExcMessage("Error: missing compress() call."));
1648 }
1649
1650
1651
1652 inline void
1653 MatrixBase::prepare_add()
1654 {
1655 prepare_action(VectorOperation::add);
1656 }
1657
1658
1659
1660 inline void
1661 MatrixBase::prepare_set()
1662 {
1663 prepare_action(VectorOperation::insert);
1664 }
1665
1666 inline MPI_Comm
1667 MatrixBase::get_mpi_communicator() const
1668 {
1669 return PetscObjectComm(reinterpret_cast<PetscObject>(matrix));
1670 }
1671
1672# endif // DOXYGEN
1673} // namespace PETScWrappers
1674
1675
1677
1678
1679#endif // DEAL_II_WITH_PETSC
1680
1681#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
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 clear_rows(const ArrayView< const size_type > &rows, const PetscScalar new_diag_value=0)
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:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
Point< 2 > second
Definition grid_out.cc:4614
Point< 2 > first
Definition grid_out.cc:4613
__global__ void set(Number *val, const Number s, const size_type N)
#define DeclException0(Exception0)
Definition exceptions.h:471
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:539
#define AssertIndexRange(index, range)
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:494
static ::ExceptionBase & ExcSourceEqualsDestination()
static ::ExceptionBase & ExcWrongMode(int arg1, int arg2)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition exceptions.h:562
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
const Iterator const std_cxx20::type_identity_t< Iterator > & end
Definition parallel.h:610
const Iterator & begin
Definition parallel.h:609
unsigned int global_dof_index
Definition types.h:81