]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Allow some matrix-vector ops to run on distributed objects.
authoryoung <young@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 8 Aug 2012 13:25:52 +0000 (13:25 +0000)
committeryoung <young@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 8 Aug 2012 13:25:52 +0000 (13:25 +0000)
git-svn-id: https://svn.dealii.org/trunk@25771 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/lac/petsc_matrix_base.h
deal.II/include/deal.II/lac/petsc_parallel_sparse_matrix.h
deal.II/include/deal.II/lac/petsc_sparse_matrix.h
deal.II/source/lac/petsc_matrix_base.cc
deal.II/source/lac/petsc_parallel_sparse_matrix.cc
deal.II/source/lac/petsc_sparse_matrix.cc

index e35ecd55bfb26137d2b23477c3deb00d77b0b912..727eee8567714bf134be7d103c10e6c6bd9fa7ed 100644 (file)
@@ -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
index 6978daf8eda5b880820762bbd7fa12027366f5c8..dae469c5286125ee57e78772945bda72c537abbd 100644 (file)
@@ -18,6 +18,7 @@
 
 #  include <deal.II/lac/exceptions.h>
 #  include <deal.II/lac/petsc_matrix_base.h>
+#  include <deal.II/lac/petsc_parallel_vector.h>
 #  include <vector>
 
 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:
 
                                          /**
index b3f21fe83c9821bd3e6854af8231d049b22922ee..aebee7087b0a7e6ed2a372036da01b4ecea4d3ba 100644 (file)
@@ -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:
 
                                        /**
index d074a21ce8106e0cff3b96b0379b53a2bc749cd8..e0e3220c2f63bb435c1ff6eebc89cb58337e58dc 100644 (file)
@@ -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,
index ef7b6c3450497a683d405d9c0a05143fcc0772db..dddd2f579a7ac45e20200ebfd31905fbe70e892d 100644 (file)
@@ -613,8 +613,6 @@ namespace PETScWrappers
         }
     }
 
-
-
     // explicit instantiations
     //
     template
@@ -679,6 +677,25 @@ namespace PETScWrappers
                              const std::vector<unsigned int> &,
                              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;
+    }
+    
   }
 }
 
index 6e15c4d5bb89663cb0e50edc897fecaa4a4d9e36..c69443649718ae4b7d16ee892eb1f31789513bf0 100644 (file)
@@ -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;
+  }
 }
 
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.