From 4699d7ffc30e85ff1b6c62852b7918905dd69532 Mon Sep 17 00:00:00 2001 From: young Date: Wed, 8 Aug 2012 13:25:52 +0000 Subject: [PATCH] Allow some matrix-vector ops to run on distributed objects. git-svn-id: https://svn.dealii.org/trunk@25771 0785d39b-7218-0410-832d-ea1e28bc413d --- .../include/deal.II/lac/petsc_matrix_base.h | 62 ------------------- .../lac/petsc_parallel_sparse_matrix.h | 48 ++++++++++++++ .../include/deal.II/lac/petsc_sparse_matrix.h | 44 +++++++++++++ deal.II/source/lac/petsc_matrix_base.cc | 22 ------- .../lac/petsc_parallel_sparse_matrix.cc | 21 ++++++- deal.II/source/lac/petsc_sparse_matrix.cc | 18 ++++++ 6 files changed, 129 insertions(+), 86 deletions(-) diff --git a/deal.II/include/deal.II/lac/petsc_matrix_base.h b/deal.II/include/deal.II/lac/petsc_matrix_base.h index e35ecd55bf..727eee8567 100644 --- a/deal.II/include/deal.II/lac/petsc_matrix_base.h +++ b/deal.II/include/deal.II/lac/petsc_matrix_base.h @@ -952,68 +952,6 @@ namespace PETScWrappers void Tvmult_add (VectorBase &dst, const VectorBase &src) 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. - * - * Note that if the current object - * represents a parallel distributed - * matrix (of type - * PETScWrappers::MPI::SparseMatrix), - * then the given vector has to be - * a distributed vector as - * well. Conversely, if the matrix is - * not distributed, then neither - * may the vector be. - */ - 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. - * - * Note that if the current object - * represents a parallel distributed - * matrix (of type - * PETScWrappers::MPI::SparseMatrix), - * then both vectors have to be - * distributed vectors as - * well. Conversely, if the matrix is - * not distributed, then neither of the - * vectors may be. - */ - PetscScalar matrix_scalar_product (const VectorBase &u, - const VectorBase &v) const; /** * Compute the residual of an diff --git a/deal.II/include/deal.II/lac/petsc_parallel_sparse_matrix.h b/deal.II/include/deal.II/lac/petsc_parallel_sparse_matrix.h index 6978daf8ed..dae469c528 100644 --- a/deal.II/include/deal.II/lac/petsc_parallel_sparse_matrix.h +++ b/deal.II/include/deal.II/lac/petsc_parallel_sparse_matrix.h @@ -18,6 +18,7 @@ # include # include +# include # include DEAL_II_NAMESPACE_OPEN @@ -32,6 +33,8 @@ namespace PETScWrappers namespace MPI { + + /** * Implementation of a parallel sparse matrix class based on PETSC, with rows * of the matrix distributed across an MPI network. All the functionality is @@ -109,6 +112,7 @@ namespace PETScWrappers class SparseMatrix : public MatrixBase { public: + /** * A structure that describes some of * the traits of this class in terms @@ -381,6 +385,50 @@ namespace PETScWrappers << "The number of local rows " << arg1 << " must be larger than the total number of rows " << arg2); //@} + + /** + * 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 Vector &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 Vector &u, + const Vector &v) const; + private: /** diff --git a/deal.II/include/deal.II/lac/petsc_sparse_matrix.h b/deal.II/include/deal.II/lac/petsc_sparse_matrix.h index b3f21fe83c..aebee7087b 100644 --- a/deal.II/include/deal.II/lac/petsc_sparse_matrix.h +++ b/deal.II/include/deal.II/lac/petsc_sparse_matrix.h @@ -44,6 +44,7 @@ namespace PETScWrappers class SparseMatrix : public MatrixBase { public: + /** * A structure that describes some of * the traits of this class in terms of @@ -275,6 +276,49 @@ 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; + private: /** diff --git a/deal.II/source/lac/petsc_matrix_base.cc b/deal.II/source/lac/petsc_matrix_base.cc index d074a21ce8..e0e3220c2f 100644 --- a/deal.II/source/lac/petsc_matrix_base.cc +++ b/deal.II/source/lac/petsc_matrix_base.cc @@ -549,28 +549,6 @@ namespace PETScWrappers } - - PetscScalar - MatrixBase::matrix_norm_square (const VectorBase &v) const - { - Vector tmp(v.size()); - vmult (tmp, v); - return tmp*v; - } - - - - PetscScalar - MatrixBase::matrix_scalar_product (const VectorBase &u, - const VectorBase &v) const - { - Vector tmp(v.size()); - vmult (tmp, v); - return u*tmp; - } - - - PetscScalar MatrixBase::residual (VectorBase &dst, const VectorBase &x, diff --git a/deal.II/source/lac/petsc_parallel_sparse_matrix.cc b/deal.II/source/lac/petsc_parallel_sparse_matrix.cc index ef7b6c3450..dddd2f579a 100644 --- a/deal.II/source/lac/petsc_parallel_sparse_matrix.cc +++ b/deal.II/source/lac/petsc_parallel_sparse_matrix.cc @@ -613,8 +613,6 @@ namespace PETScWrappers } } - - // explicit instantiations // template @@ -679,6 +677,25 @@ namespace PETScWrappers const std::vector &, const unsigned int , const bool); + + + PetscScalar + SparseMatrix::matrix_norm_square (const Vector &v) const + { + Vector tmp (v); + vmult (tmp, v); + return tmp*v; + } + + PetscScalar + SparseMatrix::matrix_scalar_product (const Vector &u, + const Vector &v) const + { + Vector tmp (v); + vmult (tmp, v); + return u*tmp; + } + } } diff --git a/deal.II/source/lac/petsc_sparse_matrix.cc b/deal.II/source/lac/petsc_sparse_matrix.cc index 6e15c4d5bb..c694436497 100644 --- a/deal.II/source/lac/petsc_sparse_matrix.cc +++ b/deal.II/source/lac/petsc_sparse_matrix.cc @@ -25,6 +25,7 @@ DEAL_II_NAMESPACE_OPEN namespace PETScWrappers { + SparseMatrix::SparseMatrix () { const int m=0, n=0, n_nonzero_per_row=0; @@ -351,6 +352,23 @@ namespace PETScWrappers template void SparseMatrix::do_reinit (const CompressedSimpleSparsityPattern &, 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