|
Reference documentation for deal.II version 9.2.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\}}\)
Go to the documentation of this file.
16 #ifndef dealii_petsc_matrix_base_h
17 # define dealii_petsc_matrix_base_h
22 # ifdef DEAL_II_WITH_PETSC
32 # include <petscmat.h>
42 template <
typename Matrix>
52 namespace MatrixIterators
125 <<
"You tried to access row " << arg1
126 <<
" of a distributed matrix, but only rows " << arg2
127 <<
" through " << arg3
128 <<
" are stored locally and can be accessed.");
239 <<
"Attempt to access element " << arg2 <<
" of row "
240 << arg1 <<
" which doesn't have that many elements.");
377 set(
const std::vector<size_type> & indices,
379 const bool elide_zero_values =
false);
387 set(
const std::vector<size_type> & row_indices,
388 const std::vector<size_type> & col_indices,
390 const bool elide_zero_values =
false);
408 const std::vector<size_type> & col_indices,
409 const std::vector<PetscScalar> &values,
410 const bool elide_zero_values =
false);
430 const PetscScalar *values,
431 const bool elide_zero_values =
false);
465 add(
const std::vector<size_type> & indices,
467 const bool elide_zero_values =
true);
475 add(
const std::vector<size_type> & row_indices,
476 const std::vector<size_type> & col_indices,
478 const bool elide_zero_values =
true);
496 const std::vector<size_type> & col_indices,
497 const std::vector<PetscScalar> &values,
498 const bool elide_zero_values =
true);
518 const PetscScalar *values,
519 const bool elide_zero_values =
true,
520 const bool col_indices_are_sorted =
false);
550 clear_rows(
const std::vector<size_type> &rows,
551 const PetscScalar new_diag_value = 0);
634 std::pair<size_type, size_type>
887 operator Mat()
const;
925 write_ascii(
const PetscViewerFormat format = PETSC_VIEWER_DEFAULT);
935 print(std::ostream &out,
const bool alternative_output =
false)
const;
947 "You are attempting an operation on two matrices that "
948 "are the same object, but the operation requires that the "
949 "two objects are in fact different.");
957 <<
"You tried to do a "
958 << (arg1 == 1 ?
"'set'" : (arg1 == 2 ?
"'add'" :
"???"))
959 <<
" operation but the matrix is currently in "
960 << (arg2 == 1 ?
"'set'" : (arg2 == 2 ?
"'add'" :
"???"))
961 <<
" mode. You first have to call 'compress()'.");
1078 friend class ::BlockMatrixBase;
1087 namespace MatrixIterators
1096 visit_present_row();
1113 return (*colnum_cache)[a_index];
1129 return (*value_cache)[a_index];
1141 inline const_iterator &
1150 if (accessor.a_index >= accessor.colnum_cache->size())
1152 accessor.a_index = 0;
1155 while ((accessor.a_row < accessor.matrix->m()) &&
1156 (accessor.a_row < accessor.matrix->local_range().second) &&
1157 (accessor.matrix->row_length(accessor.a_row) == 0))
1160 accessor.visit_present_row();
1166 inline const_iterator
1169 const const_iterator old_state = *
this;
1190 return (accessor.a_row == other.accessor.a_row &&
1191 accessor.a_index == other.accessor.a_index);
1198 return !(*
this == other);
1205 return (accessor.row() < other.accessor.row() ||
1206 (accessor.row() == other.accessor.row() &&
1207 accessor.index() < other.accessor.index()));
1234 const bool elide_zero_values)
1236 Assert(indices.size() == values.
m(),
1240 for (
size_type i = 0; i < indices.size(); ++i)
1252 const std::vector<size_type> & col_indices,
1254 const bool elide_zero_values)
1256 Assert(row_indices.size() == values.
m(),
1258 Assert(col_indices.size() == values.
n(),
1261 for (
size_type i = 0; i < row_indices.size(); ++i)
1273 const std::vector<size_type> & col_indices,
1274 const std::vector<PetscScalar> &values,
1275 const bool elide_zero_values)
1277 Assert(col_indices.size() == values.size(),
1293 const PetscScalar *values,
1294 const bool elide_zero_values)
1298 const PetscInt petsc_i = row;
1299 PetscInt
const *col_index_ptr;
1301 PetscScalar
const *col_value_ptr;
1305 if (elide_zero_values ==
false)
1307 col_index_ptr =
reinterpret_cast<const PetscInt *
>(col_indices);
1308 col_value_ptr = values;
1315 if (column_indices.size() < n_cols)
1317 column_indices.resize(n_cols);
1318 column_values.resize(n_cols);
1324 const PetscScalar
value = values[j];
1326 if (
value != PetscScalar())
1328 column_indices[n_columns] = col_indices[j];
1329 column_values[n_columns] =
value;
1335 col_index_ptr = column_indices.data();
1336 col_value_ptr = column_values.data();
1339 const PetscErrorCode ierr = MatSetValues(
matrix,
1356 if (
value == PetscScalar())
1367 add(i, 1, &j, &
value,
false);
1373 MatrixBase::add(
const std::vector<size_type> & indices,
1375 const bool elide_zero_values)
1377 Assert(indices.size() == values.
m(),
1381 for (
size_type i = 0; i < indices.size(); ++i)
1392 MatrixBase::add(
const std::vector<size_type> & row_indices,
1393 const std::vector<size_type> & col_indices,
1395 const bool elide_zero_values)
1397 Assert(row_indices.size() == values.
m(),
1399 Assert(col_indices.size() == values.
n(),
1402 for (
size_type i = 0; i < row_indices.size(); ++i)
1414 const std::vector<size_type> & col_indices,
1415 const std::vector<PetscScalar> &values,
1416 const bool elide_zero_values)
1418 Assert(col_indices.size() == values.size(),
1434 const PetscScalar *values,
1435 const bool elide_zero_values,
1438 (void)elide_zero_values;
1442 const PetscInt petsc_i = row;
1443 PetscInt
const *col_index_ptr;
1445 PetscScalar
const *col_value_ptr;
1449 if (elide_zero_values ==
false)
1451 col_index_ptr =
reinterpret_cast<const PetscInt *
>(col_indices);
1452 col_value_ptr = values;
1459 if (column_indices.size() < n_cols)
1461 column_indices.resize(n_cols);
1462 column_values.resize(n_cols);
1468 const PetscScalar
value = values[j];
1470 if (
value != PetscScalar())
1472 column_indices[n_columns] = col_indices[j];
1473 column_values[n_columns] =
value;
1479 col_index_ptr = column_indices.data();
1480 col_value_ptr = column_values.data();
1483 const PetscErrorCode ierr = MatSetValues(
1484 matrix, 1, &petsc_i, n_columns, col_index_ptr, col_value_ptr, ADD_VALUES);
1498 inline MatrixBase::const_iterator
1502 (in_local_range(0) && in_local_range(m() - 1)),
1504 "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."));
1509 while ((first_nonempty_row < m()) && (row_length(first_nonempty_row) == 0))
1510 ++first_nonempty_row;
1512 return const_iterator(
this, first_nonempty_row, 0);
1516 inline MatrixBase::const_iterator
1520 (in_local_range(0) && in_local_range(m() - 1)),
1522 "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."));
1524 return const_iterator(
this, m(), 0);
1528 inline MatrixBase::const_iterator
1531 Assert(in_local_range(r),
1534 if (row_length(r) > 0)
1535 return const_iterator(
this, r, 0);
1541 inline MatrixBase::const_iterator
1544 Assert(in_local_range(r),
1558 if (i == local_range().second || (row_length(i) > 0))
1559 return const_iterator(
this, i, 0);
1565 return {
this, m(), 0};
1571 MatrixBase::in_local_range(
const size_type index)
const
1575 const PetscErrorCode ierr =
1576 MatGetOwnershipRange(
static_cast<const Mat &
>(
matrix), &
begin, &
end);
1589 last_action = new_action;
1591 Assert(last_action == new_action, ExcWrongMode(last_action, new_action));
1597 MatrixBase::assert_is_compressed()
1602 ExcMessage(
"Error: missing compress() call."));
1608 MatrixBase::prepare_add()
1616 MatrixBase::prepare_set()
1628 # endif // DEAL_II_WITH_PETSC
std::size_t memory_consumption() const
std::vector< PetscInt > column_indices
#define DeclExceptionMsg(Exception, defaulttext)
bool operator!=(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
types::global_dof_index size_type
std::shared_ptr< const std::vector< size_type > > colnum_cache
SynchronousIterators< Iterators > operator++(SynchronousIterators< Iterators > &a)
bool operator<(const SynchronousIterators< Iterators > &a, const SynchronousIterators< Iterators > &b)
void Tmmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const
void vmult_add(VectorBase &dst, const VectorBase &src) const
bool in_local_range(const size_type index) const
const_iterator(const MatrixBase *matrix, const size_type row, const size_type index)
static ::ExceptionBase & ExcBeyondEndOfMatrix()
bool operator==(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
MatrixTableIterators::Accessor< TransposeTable< T >, Constness, MatrixTableIterators::Storage::column_major > Accessor
void mmult(MatrixBase &C, const MatrixBase &B, const VectorBase &V) const
PetscScalar matrix_scalar_product(const VectorBase &u, const VectorBase &v) const
const Accessor & operator*() const
PetscScalar matrix_norm_square(const VectorBase &v) const
size_type n_nonzero_elements() const
types::global_dof_index size_type
static ::ExceptionBase & ExcInvalidIndexWithinRow(int arg1, int arg2)
size_type local_size() const
PetscReal l1_norm() const
void set(const size_type i, const size_type j, const PetscScalar value)
#define AssertIndexRange(index, range)
SymmetricTensor< 2, dim, Number > C(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
const_iterator & operator++()
void print(std::ostream &out, const bool alternative_output=false) const
bool operator!=(const const_iterator &) const
PetscBool is_symmetric(const double tolerance=1.e-12)
void prepare_action(const VectorOperation::values new_action)
std::shared_ptr< const std::vector< PetscScalar > > value_cache
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
Accessor(const MatrixBase *matrix, const size_type row, const size_type index)
MatrixBase & operator*=(const PetscScalar factor)
__global__ void set(Number *val, const Number s, const size_type N)
void vmult(VectorBase &dst, const VectorBase &src) const
std::vector< PetscScalar > column_values
PetscReal frobenius_norm() const
static ::ExceptionBase & ExcNotQuadratic()
PetscScalar diag_element(const size_type i) const
size_type row_length(const size_type row) const
PetscScalar trace() const
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
PetscScalar value() const
#define AssertIsFinite(number)
void Tvmult(VectorBase &dst, const VectorBase &src) const
void write_ascii(const PetscViewerFormat format=PETSC_VIEWER_DEFAULT)
static ::ExceptionBase & ExcMessage(std::string arg1)
unsigned int global_dof_index
VectorType::value_type * begin(VectorType &V)
bool operator==(const const_iterator &) const
#define DEAL_II_NAMESPACE_OPEN
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
std::pair< size_type, size_type > local_range() const
virtual const MPI_Comm & get_mpi_communicator() const =0
@ matrix
Contents is actually a matrix.
void clear_rows(const std::vector< size_type > &rows, const PetscScalar new_diag_value=0)
void add(const size_type i, const size_type j, const PetscScalar value)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
MatrixBase & operator/=(const PetscScalar factor)
PetscBool is_hermitian(const double tolerance=1.e-12)
VectorType::value_type * end(VectorType &V)
#define DEAL_II_DEPRECATED
void compress(const VectorOperation::values operation)
void assert_is_compressed()
types::global_dof_index size_type
#define Assert(cond, exc)
void clear_row(const size_type row, const PetscScalar new_diag_value=0)
#define DeclException0(Exception0)
PetscScalar residual(VectorBase &dst, const VectorBase &x, const VectorBase &b) const
void Tvmult_add(VectorBase &dst, const VectorBase &src) const
const_iterator end() const
static ::ExceptionBase & ExcSourceEqualsDestination()
const_iterator begin() const
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
bool operator<(const const_iterator &) const
MatrixBase & operator=(const MatrixBase &)=delete
PetscScalar el(const size_type i, const size_type j) const
static ::ExceptionBase & ExcWrongMode(int arg1, int arg2)
virtual ~MatrixBase() override
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcAccessToNonlocalRow(int arg1, int arg2, int arg3)
VectorOperation::values last_action
PetscReal linfty_norm() const
#define DeclException2(Exception2, type1, type2, outsequence)
#define AssertThrow(cond, exc)
static ::ExceptionBase & ExcIteratorPastEnd()
std::enable_if< std::is_floating_point< T >::value &&std::is_floating_point< U >::value, typename ProductType< std::complex< T >, std::complex< U > >::type >::type operator*(const std::complex< T > &left, const std::complex< U > &right)
PetscScalar operator()(const size_type i, const size_type j) const
const Accessor * operator->() const
static ::ExceptionBase & ExcBeyondEndOfMatrix()