17#ifdef DEAL_II_WITH_PETSC
29 namespace MatrixIterators
33 MatrixBase::const_iterator::Accessor::visit_present_row()
38 if (
matrix->in_local_range(this->a_row) ==
false)
48 const PetscInt *colnums;
52 MatGetRow(*matrix, this->a_row, &ncols, &colnums, &values);
64 for (PetscInt j = 0; j < ncols; ++j)
66 const auto column =
static_cast<PetscInt
>(colnums[j]);
71 std::make_shared<std::vector<size_type>>(colnums, colnums + ncols);
73 std::make_shared<std::vector<PetscScalar>>(
values,
values + ncols);
76 ierr = MatRestoreRow(*matrix, this->a_row, &ncols, &colnums, &values);
94 const PetscErrorCode ierr =
95 PetscObjectReference(
reinterpret_cast<PetscObject
>(
matrix));
104 PetscErrorCode ierr =
105 PetscObjectReference(
reinterpret_cast<PetscObject
>(A));
107 ierr = MatDestroy(&
matrix);
114 PetscErrorCode ierr = MatDestroy(&
matrix);
124 const PetscErrorCode ierr = MatDestroy(&
matrix);
130 const int m = 0,
n = 0, n_nonzero_per_row = 0;
131 const PetscErrorCode ierr = MatCreateSeqAIJ(
132 PETSC_COMM_SELF,
m,
n, n_nonzero_per_row,
nullptr, &
matrix);
146 const PetscErrorCode ierr = MatZeroEntries(
matrix);
164 const PetscScalar new_diag_value)
171 for (
const auto &row : rows)
174 const std::vector<PetscInt> petsc_rows(rows.
begin(), rows.
end());
189 ierr = MatZeroRowsIS(
matrix, index_set, new_diag_value,
nullptr,
nullptr);
191 ierr = ISDestroy(&index_set);
197 const PetscScalar new_diag_value)
204 for (
const auto &row : rows)
207 const std::vector<PetscInt> petsc_rows(rows.begin(), rows.end());
223 MatZeroRowsColumnsIS(
matrix, index_set, new_diag_value,
nullptr,
nullptr);
225 ierr = ISDestroy(&index_set);
234 const auto petsc_i =
static_cast<PetscInt
>(i);
236 const auto petsc_j =
static_cast<PetscInt
>(j);
241 const PetscErrorCode ierr =
242 MatGetValues(
matrix, 1, &petsc_i, 1, &petsc_j, &value);
272 int all_int_last_action;
274 const int ierr = MPI_Allreduce(&my_int_last_action,
275 &all_int_last_action,
285 "Error: not all processors agree on the last "
286 "VectorOperation before this compress() call."));
293 "Missing compress() or calling with wrong VectorOperation argument."));
296 PetscErrorCode ierr = MatAssemblyBegin(
matrix, MAT_FINAL_ASSEMBLY);
299 ierr = MatAssemblyEnd(
matrix, MAT_FINAL_ASSEMBLY);
310 PetscInt n_rows, n_cols;
312 const PetscErrorCode ierr = MatGetSize(
matrix, &n_rows, &n_cols);
323 PetscInt n_rows, n_cols;
325 const PetscErrorCode ierr = MatGetSize(
matrix, &n_rows, &n_cols);
338 const PetscErrorCode ierr = MatGetLocalSize(
matrix, &n_rows,
nullptr);
346 std::pair<MatrixBase::size_type, MatrixBase::size_type>
351 const PetscErrorCode ierr =
352 MatGetOwnershipRange(
static_cast<const Mat &
>(
matrix), &
begin, &
end);
365 const PetscErrorCode ierr = MatGetLocalSize(
matrix,
nullptr, &n_cols);
373 std::pair<MatrixBase::size_type, MatrixBase::size_type>
378 const PetscErrorCode ierr =
379 MatGetOwnershipRangeColumn(
static_cast<const Mat &
>(
matrix),
393 const PetscErrorCode ierr = MatGetInfo(
matrix, MAT_GLOBAL_SUM, &mat_info);
398 return static_cast<std::uint64_t
>(mat_info.nz_used);
418 const PetscInt *colnums;
419 const PetscScalar *values;
423 PetscErrorCode ierr = MatGetRow(*
this, row, &ncols, &colnums, &values);
431 const PetscInt ncols_saved = ncols;
432 ierr = MatRestoreRow(*
this, row, &ncols, &colnums, &values);
444 const PetscErrorCode ierr = MatNorm(
matrix, NORM_1, &result);
457 const PetscErrorCode ierr = MatNorm(
matrix, NORM_INFINITY, &result);
470 const PetscErrorCode ierr = MatNorm(
matrix, NORM_FROBENIUS, &result);
506 const PetscErrorCode ierr = MatGetTrace(
matrix, &result);
517 const PetscErrorCode ierr = MatScale(
matrix, a);
528 const PetscScalar factor = 1. / a;
529 const PetscErrorCode ierr = MatScale(
matrix, factor);
540 const PetscErrorCode ierr =
541 MatAXPY(
matrix, factor, other, DIFFERENT_NONZERO_PATTERN);
553 const PetscErrorCode ierr = MatMult(
matrix, src, dst);
564 const PetscErrorCode ierr = MatMultTranspose(
matrix, src, dst);
575 const PetscErrorCode ierr = MatMultAdd(
matrix, src, dst, dst);
586 const PetscErrorCode ierr = MatMultTransposeAdd(
matrix, src, dst, dst);
598 const bool transpose_left)
600 const bool use_vector = (V.size() == inputright.
m() ? true :
false);
601 if (transpose_left ==
false)
603 Assert(inputleft.
n() == inputright.
m(),
608 Assert(inputleft.
m() == inputright.
m(),
620 ierr = MatTransposeMatMult(inputleft,
629 ierr = MatMatMult(inputleft,
640 ierr = MatDuplicate(inputleft, MAT_COPY_VALUES, &tmp);
644# if DEAL_II_PETSC_VERSION_LT(3, 8, 0)
645 ierr = MatTranspose(tmp, MAT_REUSE_MATRIX, &tmp);
647 ierr = MatTranspose(tmp, MAT_INPLACE_MATRIX, &tmp);
651 ierr = MatDiagonalScale(tmp,
nullptr, V);
653 ierr = MatMatMult(tmp,
659 ierr = MatDestroy(&tmp);
698 MatrixBase::operator Mat()
const
698 MatrixBase::operator Mat()
const {
…}
712# if DEAL_II_PETSC_VERSION_LT(3, 8, 0)
713 const PetscErrorCode ierr = MatTranspose(
matrix, MAT_REUSE_MATRIX, &
matrix);
715 const PetscErrorCode ierr =
726 const PetscErrorCode ierr = MatIsSymmetric(
matrix, tolerance, &truth);
737 const PetscErrorCode ierr = MatIsHermitian(
matrix, tolerance, &truth);
750 PetscErrorCode ierr =
751 PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(
comm), format);
755 ierr = MatView(
matrix, PETSC_VIEWER_STDOUT_(
comm));
757 ierr = PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(
comm));
766 PetscErrorCode ierr = MatHasOperation(
matrix, MATOP_GET_ROW, &has);
772 ierr = MatConvert(
matrix, MATAIJ, MAT_INITIAL_MATRIX, &vmatrix);
776 std::pair<MatrixBase::size_type, MatrixBase::size_type> loc_range =
780 const PetscInt *colnums;
781 const PetscScalar *values;
784 for (row = loc_range.first; row < loc_range.second; ++row)
786 ierr = MatGetRow(vmatrix, row, &ncols, &colnums, &values);
789 for (PetscInt col = 0; col < ncols; ++col)
791 out <<
"(" << row <<
"," << colnums[col] <<
") " << values[col]
795 ierr = MatRestoreRow(vmatrix, row, &ncols, &colnums, &values);
800 ierr = MatDestroy(&vmatrix);
812 const PetscErrorCode ierr = MatGetInfo(
matrix, MAT_LOCAL, &info);
815 return sizeof(*this) +
static_cast<size_type>(info.memory);
void add(const size_type i, const size_type j, const PetscScalar value)
size_type row_length(const size_type row) const
std::size_t memory_consumption() const
PetscReal l1_norm() const
VectorOperation::values last_action
void vmult(VectorBase &dst, const VectorBase &src) const
MPI_Comm get_mpi_communicator() const
PetscScalar diag_element(const size_type i) const
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)
PetscReal frobenius_norm() const
virtual ~MatrixBase() override
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
PetscScalar el(const size_type i, const size_type j) const
size_type local_size() const
MatrixBase & operator=(const MatrixBase &)=delete
PetscScalar trace() const
PetscBool is_hermitian(const double tolerance=1.e-12)
PetscScalar matrix_scalar_product(const VectorBase &u, const VectorBase &v) const
void assert_is_compressed()
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 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 clear_rows(const ArrayView< const size_type > &rows, const PetscScalar new_diag_value=0)
std::uint64_t n_nonzero_elements() const
void clear_row(const size_type row, const PetscScalar new_diag_value=0)
PetscReal linfty_norm() const
real_type l2_norm() const
size_type size() const override
#define DEAL_II_NAMESPACE_OPEN
constexpr bool running_in_debug_mode()
#define DEAL_II_NAMESPACE_CLOSE
#define AssertIntegerConversion(index1, index2)
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcScalarAssignmentOnlyForZeroValue()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertThrowMPI(error_code)
#define AssertNothrow(cond, exc)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcSourceEqualsDestination()
static ::ExceptionBase & ExcInternalError()
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.
void perform_mmult(const MatrixBase &inputleft, const MatrixBase &inputright, MatrixBase &result, const VectorBase &V, const bool transpose_left)