From c90fd9ff8aa5e06abe911c9ddd453010335ddcaa Mon Sep 17 00:00:00 2001 From: David Wells Date: Sat, 29 Apr 2017 12:35:55 -0400 Subject: [PATCH] Get rid of PETScWrappers::Vector in the PETSc matrix classes. We can simply use PETScWrappers::VectorBase instead. In addition, get rid of some methods in the derived class that are identical to those in MatrixBase. --- include/deal.II/lac/petsc_sparse_matrix.h | 27 ----------------------- source/lac/petsc_matrix_base.cc | 10 ++++----- source/lac/petsc_sparse_matrix.cc | 19 +--------------- 3 files changed, 6 insertions(+), 50 deletions(-) diff --git a/include/deal.II/lac/petsc_sparse_matrix.h b/include/deal.II/lac/petsc_sparse_matrix.h index 41a1481f71..ce88ecc780 100644 --- a/include/deal.II/lac/petsc_sparse_matrix.h +++ b/include/deal.II/lac/petsc_sparse_matrix.h @@ -198,33 +198,6 @@ namespace PETScWrappers */ virtual const MPI_Comm &get_mpi_communicator () const; - /** - * Return the square of the norm of the vector $v$ with respect to the - * norm induced by this matrix, i.e. $\left(v,Mv\right)$. This is useful, - * e.g. in the finite element context, where the $L_2$ norm of a function - * equals the matrix norm with respect to the mass matrix of the vector - * representing the nodal values of the finite element function. - * - * Obviously, the matrix needs to be quadratic for this operation. - * - * The implementation of this function is not as efficient as the one in - * the @p MatrixBase class used in deal.II (i.e. the original one, not the - * PETSc wrapper class) since PETSc doesn't support this operation and - * needs a temporary vector. - */ - PetscScalar matrix_norm_square (const VectorBase &v) const; - - /** - * Compute the matrix scalar product $\left(u,Mv\right)$. - * - * The implementation of this function is not as efficient as the one in - * the @p MatrixBase class used in deal.II (i.e. the original one, not the - * PETSc wrapper class) since PETSc doesn't support this operation and - * needs a temporary vector. - */ - PetscScalar matrix_scalar_product (const VectorBase &u, - const VectorBase &v) const; - /** * Return the number of rows of this matrix. */ diff --git a/source/lac/petsc_matrix_base.cc b/source/lac/petsc_matrix_base.cc index ac5c4b77db..41dc140cd1 100644 --- a/source/lac/petsc_matrix_base.cc +++ b/source/lac/petsc_matrix_base.cc @@ -22,7 +22,7 @@ # include # include # include -# include +# include DEAL_II_NAMESPACE_OPEN @@ -384,8 +384,8 @@ namespace PETScWrappers PetscScalar MatrixBase::matrix_norm_square (const VectorBase &v) const { - Vector tmp(v.size()); - vmult (tmp, v); + VectorBase tmp(v); + vmult(tmp, v); return tmp*v; } @@ -394,8 +394,8 @@ namespace PETScWrappers MatrixBase::matrix_scalar_product (const VectorBase &u, const VectorBase &v) const { - Vector tmp(v.size()); - vmult (tmp, v); + VectorBase tmp(u); + vmult(tmp, v); return u*tmp; } diff --git a/source/lac/petsc_sparse_matrix.cc b/source/lac/petsc_sparse_matrix.cc index 913490b7e8..e71099b712 100644 --- a/source/lac/petsc_sparse_matrix.cc +++ b/source/lac/petsc_sparse_matrix.cc @@ -19,7 +19,7 @@ # include # include -# include +# include # include # include @@ -284,23 +284,6 @@ namespace PETScWrappers template void SparseMatrix::do_reinit (const DynamicSparsityPattern &, const bool); - - PetscScalar - SparseMatrix::matrix_norm_square (const VectorBase &v) const - { - Vector tmp (v.size()); - vmult (tmp, v); - return tmp*v; - } - - PetscScalar - SparseMatrix::matrix_scalar_product (const VectorBase &u, - const VectorBase &v) const - { - Vector tmp (v.size()); - vmult (tmp, v); - return u*tmp; - } } -- 2.39.5