18 #ifdef DEAL_II_WITH_PETSC
30 namespace MatrixIterators
39 if (
matrix->in_local_range(this->a_row) ==
false)
49 const PetscInt * colnums;
53 MatGetRow(*
matrix, this->a_row, &ncols, &colnums, &
values);
64 std::make_shared<std::vector<size_type>>(colnums, colnums + ncols);
66 std::make_shared<std::vector<PetscScalar>>(
values,
values + ncols);
69 ierr = MatRestoreRow(*
matrix, this->a_row, &ncols, &colnums, &
values);
87 const PetscErrorCode ierr =
88 PetscObjectReference(
reinterpret_cast<PetscObject
>(
matrix));
98 PetscObjectReference(
reinterpret_cast<PetscObject
>(
A));
120 const int m = 0,
n = 0, n_nonzero_per_row = 0;
121 const PetscErrorCode ierr = MatCreateSeqAIJ(
122 PETSC_COMM_SELF,
m,
n, n_nonzero_per_row,
nullptr, &
matrix);
136 const PetscErrorCode ierr = MatZeroEntries(
matrix);
147 std::vector<size_type> rows(1, row);
155 const PetscScalar new_diag_value)
161 const std::vector<PetscInt> petsc_rows(rows.begin(), rows.end());
174 const PetscErrorCode ierr =
175 MatZeroRowsIS(
matrix, index_set, new_diag_value,
nullptr,
nullptr);
177 ISDestroy(&index_set);
182 const PetscScalar new_diag_value)
188 const std::vector<PetscInt> petsc_rows(rows.begin(), rows.end());
201 const PetscErrorCode ierr =
202 MatZeroRowsColumnsIS(
matrix, index_set, new_diag_value,
nullptr,
nullptr);
204 ISDestroy(&index_set);
212 PetscInt petsc_i = i, petsc_j = j;
216 const PetscErrorCode ierr =
217 MatGetValues(
matrix, 1, &petsc_i, 1, &petsc_j, &value);
242 # ifdef DEAL_II_WITH_MPI
246 int all_int_last_action;
248 const int ierr = MPI_Allreduce(&my_int_last_action,
249 &all_int_last_action,
258 ExcMessage(
"Error: not all processors agree on the last "
259 "VectorOperation before this compress() call."));
267 "Missing compress() or calling with wrong VectorOperation argument."));
270 PetscErrorCode ierr = MatAssemblyBegin(
matrix, MAT_FINAL_ASSEMBLY);
273 ierr = MatAssemblyEnd(
matrix, MAT_FINAL_ASSEMBLY);
284 PetscInt n_rows, n_cols;
286 const PetscErrorCode ierr = MatGetSize(
matrix, &n_rows, &n_cols);
297 PetscInt n_rows, n_cols;
299 const PetscErrorCode ierr = MatGetSize(
matrix, &n_rows, &n_cols);
310 PetscInt n_rows, n_cols;
312 const PetscErrorCode ierr = MatGetLocalSize(
matrix, &n_rows, &n_cols);
320 std::pair<MatrixBase::size_type, MatrixBase::size_type>
325 const PetscErrorCode ierr =
326 MatGetOwnershipRange(
static_cast<const Mat &
>(
matrix), &
begin, &
end);
338 const PetscErrorCode ierr = MatGetInfo(
matrix, MAT_GLOBAL_SUM, &mat_info);
343 return static_cast<std::uint64_t
>(mat_info.nz_used);
363 const PetscInt * colnums;
364 const PetscScalar *
values;
368 PetscErrorCode ierr = MatGetRow(*
this, row, &ncols, &colnums, &
values);
376 const PetscInt ncols_saved = ncols;
377 ierr = MatRestoreRow(*
this, row, &ncols, &colnums, &
values);
389 const PetscErrorCode ierr = MatNorm(
matrix, NORM_1, &result);
402 const PetscErrorCode ierr = MatNorm(
matrix, NORM_INFINITY, &result);
415 const PetscErrorCode ierr = MatNorm(
matrix, NORM_FROBENIUS, &result);
446 const PetscErrorCode ierr = MatGetTrace(
matrix, &result);
457 const PetscErrorCode ierr = MatScale(
matrix, a);
468 const PetscScalar factor = 1. / a;
469 const PetscErrorCode ierr = MatScale(
matrix, factor);
480 const PetscErrorCode ierr =
481 MatAXPY(
matrix, factor, other, DIFFERENT_NONZERO_PATTERN);
493 const PetscErrorCode ierr = MatMult(
matrix, src, dst);
504 const PetscErrorCode ierr = MatMultTranspose(
matrix, src, dst);
515 const PetscErrorCode ierr = MatMultAdd(
matrix, src, dst, dst);
526 const PetscErrorCode ierr = MatMultTransposeAdd(
matrix, src, dst, dst);
538 const bool transpose_left)
540 const bool use_vector = (
V.size() == inputright.
m() ? true :
false);
541 if (transpose_left ==
false)
543 Assert(inputleft.
n() == inputright.
m(),
548 Assert(inputleft.
m() == inputright.
m(),
556 if (use_vector ==
false)
560 ierr = MatTransposeMatMult(inputleft,
569 ierr = MatMatMult(inputleft,
580 ierr = MatDuplicate(inputleft, MAT_COPY_VALUES, &tmp);
584 # if DEAL_II_PETSC_VERSION_LT(3, 8, 0)
585 ierr = MatTranspose(tmp, MAT_REUSE_MATRIX, &tmp);
587 ierr = MatTranspose(tmp, MAT_INPLACE_MATRIX, &tmp);
591 ierr = MatDiagonalScale(tmp,
nullptr,
V);
593 ierr = MatMatMult(tmp,
638 MatrixBase::operator Mat()
const
652 # if DEAL_II_PETSC_VERSION_LT(3, 8, 0)
653 const PetscErrorCode ierr = MatTranspose(
matrix, MAT_REUSE_MATRIX, &
matrix);
655 const PetscErrorCode ierr =
666 const PetscErrorCode ierr = MatIsSymmetric(
matrix, tolerance, &truth);
677 const PetscErrorCode ierr = MatIsHermitian(
matrix, tolerance, &truth);
687 MPI_Comm
comm = PetscObjectComm(
reinterpret_cast<PetscObject
>(
matrix));
690 PetscErrorCode ierr =
691 PetscViewerSetFormat(PETSC_VIEWER_STDOUT_(
comm), format);
695 ierr = MatView(
matrix, PETSC_VIEWER_STDOUT_(
comm));
704 PetscErrorCode ierr = MatHasOperation(
matrix, MATOP_GET_ROW, &has);
710 ierr = MatConvert(
matrix, MATAIJ, MAT_INITIAL_MATRIX, &vmatrix);
714 std::pair<MatrixBase::size_type, MatrixBase::size_type> loc_range =
718 const PetscInt * colnums;
719 const PetscScalar *
values;
722 for (row = loc_range.first; row < loc_range.second; ++row)
724 ierr = MatGetRow(vmatrix, row, &ncols, &colnums, &
values);
727 for (PetscInt col = 0; col < ncols; ++col)
729 out <<
"(" << row <<
"," << colnums[col] <<
") " <<
values[col]
733 ierr = MatRestoreRow(vmatrix, row, &ncols, &colnums, &
values);
750 const PetscErrorCode ierr = MatGetInfo(
matrix, MAT_LOCAL, &info);
753 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
void clear_rows(const std::vector< size_type > &rows, const PetscScalar new_diag_value=0)
PetscScalar diag_element(const size_type i) 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)
types::global_dof_index size_type
PetscReal frobenius_norm() const
MatrixBase & operator=(const MatrixBase &)=delete
virtual ~MatrixBase() override
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
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)
const MPI_Comm & get_mpi_communicator() const
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
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
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcScalarAssignmentOnlyForZeroValue()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
#define AssertThrowMPI(error_code)
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcSourceEqualsDestination()
static ::ExceptionBase & ExcMessage(std::string arg1)
static ::ExceptionBase & ExcNotQuadratic()
#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)
PetscErrorCode destroy_matrix(Mat &matrix)
SymmetricTensor< 2, dim, Number > C(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)