]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Name these classes the same as those in PETScWrappers, but move them into a namespace...
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 23 Mar 2004 19:58:56 +0000 (19:58 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 23 Mar 2004 19:58:56 +0000 (19:58 +0000)
git-svn-id: https://svn.dealii.org/trunk@8852 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/petsc_parallel_vector.h
deal.II/lac/source/petsc_parallel_vector.cc

index 9faf31ef19a600c496e06af0f119c0330ca091b2..cd72d93d65c7423e1a9eb4af0a5ce12904c12156 100644 (file)
 
 namespace PETScWrappers
 {
+/**
+ * Namespace for PETSc classes that work in parallel over MPI, such as
+ * distributed vectors and matrices.
+ *
+ * @author Wolfgang Bangerth, 2004
+ */
+  namespace MPI
+  {
+    
 /**
  * Implementation of a parallel vector class based on PETSC and using MPI
  * communication to synchronise distributed operations. All the functionality
@@ -38,132 +47,134 @@ namespace PETScWrappers
  *
  * @author Wolfgang Bangerth, 2004
  */
-  class ParallelVector : public VectorBase
-  {
-    public:
-                                       /**
-                                        * Default constructor. Initialize the
-                                        * vector as empty.
-                                        */
-      ParallelVector ();
+    class Vector : public VectorBase
+    {
+      public:
+                                         /**
+                                          * Default constructor. Initialize the
+                                          * vector as empty.
+                                          */
+        Vector ();
       
-                                       /**
-                                        * Constructor. Set dimension to
-                                        * @p{n} and initialize all
-                                        * elements with zero.
-                                        *
-                                        * The constructor is made explicit to
-                                        * avoid accidents like this:
-                                        * @p{v=0;}. Presumably, the user wants
-                                        * to set every element of the vector to
-                                        * zero, but instead, what happens is
-                                        * this call: @p{v=Vector<number>(0);},
-                                        * i.e. the vector is replaced by one of
-                                        * length zero.
-                                        */
-      explicit ParallelVector (const unsigned int n,
-                               const unsigned int local_size,
-                               const MPI_Comm    &communicator);
+                                         /**
+                                          * Constructor. Set dimension to
+                                          * @p{n} and initialize all
+                                          * elements with zero.
+                                          *
+                                          * The constructor is made explicit to
+                                          * avoid accidents like this:
+                                          * @p{v=0;}. Presumably, the user wants
+                                          * to set every element of the vector to
+                                          * zero, but instead, what happens is
+                                          * this call: @p{v=Vector<number>(0);},
+                                          * i.e. the vector is replaced by one of
+                                          * length zero.
+                                          */
+        explicit Vector (const unsigned int n,
+                         const unsigned int local_size,
+                         const MPI_Comm    &communicator);
     
-                                       /**
-                                        * Copy-constructor from deal.II
-                                        * vectors. Sets the dimension to that
-                                        * of the given vector, and copies all
-                                        * elements.
-                                        */
-      template <typename Number>
-      explicit ParallelVector (const ::Vector<Number> &v,
-                               const unsigned int      local_size,
-                               const MPI_Comm         &communicator);
-
-                                       /**
-                                        * Copy-constructor the
-                                        * values from a PETSc wrapper vector
-                                        * class.
-                                        */
-      explicit ParallelVector (const VectorBase &v,
-                               const unsigned int local_size,
-                               const MPI_Comm    &communicator);
-
-                                       /**
-                                        * Set all components of the vector to
-                                        * the given number @p{s}. Simply pass
-                                        * this down to the base class, but we
-                                        * still need to declare this function
-                                        * to make the example given in the
-                                        * discussion about making the
-                                        * constructor explicit work.
-                                        */
-      ParallelVector & operator = (const PetscScalar s);
-
-                                       /**
-                                        * Copy the values of a deal.II vector
-                                        * (as opposed to those of the PETSc
-                                        * vector wrapper class) into this
-                                        * object.
-                                        */
-      template <typename number>
-      ParallelVector & operator = (const ::Vector<number> &v);
+                                         /**
+                                          * Copy-constructor from deal.II
+                                          * vectors. Sets the dimension to that
+                                          * of the given vector, and copies all
+                                          * elements.
+                                          */
+        template <typename Number>
+        explicit Vector (const ::Vector<Number> &v,
+                         const unsigned int      local_size,
+                         const MPI_Comm         &communicator);
+
+                                         /**
+                                          * Copy-constructor the
+                                          * values from a PETSc wrapper vector
+                                          * class.
+                                          */
+        explicit Vector (const VectorBase &v,
+                         const unsigned int local_size,
+                         const MPI_Comm    &communicator);
+
+                                         /**
+                                          * Set all components of the vector to
+                                          * the given number @p{s}. Simply pass
+                                          * this down to the base class, but we
+                                          * still need to declare this function
+                                          * to make the example given in the
+                                          * discussion about making the
+                                          * constructor explicit work.
+                                          */
+        Vector & operator = (const PetscScalar s);
+
+                                         /**
+                                          * Copy the values of a deal.II vector
+                                          * (as opposed to those of the PETSc
+                                          * vector wrapper class) into this
+                                          * object.
+                                          */
+        template <typename number>
+        Vector & operator = (const ::Vector<number> &v);
       
-    protected:
-                                       /**
-                                        * Create a vector of length @p{n}. For
-                                        * this class, we create a parallel
-                                        * vector. @arg n denotes the total
-                                        * size of the vector to be
-                                        * created. @arg local_size denotes how
-                                        * many of these elements shall be
-                                        * stored locally. The last argument is
-                                        * ignored for sequential vectors.
-                                        */
-      virtual void create_vector (const unsigned int n,
-                                  const unsigned int local_size);
-
-    private:
-                                       /**
-                                        * Copy of the communicator object to
-                                        * be used for this parallel vector.
-                                        */
-      MPI_Comm communicator;
-  };
+      protected:
+                                         /**
+                                          * Create a vector of length @p{n}. For
+                                          * this class, we create a parallel
+                                          * vector. @arg n denotes the total
+                                          * size of the vector to be
+                                          * created. @arg local_size denotes how
+                                          * many of these elements shall be
+                                          * stored locally. The last argument is
+                                          * ignored for sequential vectors.
+                                          */
+        virtual void create_vector (const unsigned int n,
+                                    const unsigned int local_size);
+
+      private:
+                                         /**
+                                          * Copy of the communicator object to
+                                          * be used for this parallel vector.
+                                          */
+        MPI_Comm communicator;
+    };
 
 
 
 // ------------------ template and inline functions -------------
 
 
-  template <typename number>
-  ParallelVector::ParallelVector (const ::Vector<number> &v,
-                                  const unsigned int      local_size,
-                                  const MPI_Comm         &communicator)
-                  :
-                  communicator (communicator)
-  {
-    ParallelVector::create_vector (v.size(), local_size);
+    template <typename number>
+    Vector::Vector (const ::Vector<number> &v,
+                    const unsigned int      local_size,
+                    const MPI_Comm         &communicator)
+                    :
+                    communicator (communicator)
+    {
+      Vector::create_vector (v.size(), local_size);
 
-    VectorBase::operator = (v);
-  }
+      VectorBase::operator = (v);
+    }
 
   
   
-  inline
-  ParallelVector &
-  ParallelVector::operator = (const PetscScalar s)
-  {
-    VectorBase::operator = (s);
-
-    return *this;
-  }
+    inline
+    Vector &
+    Vector::operator = (const PetscScalar s)
+    {
+      VectorBase::operator = (s);
+
+      return *this;
+    }
   
 
-  template <typename number>
-  inline
-  ParallelVector &
-  ParallelVector::operator = (const ::Vector<number> &v)
-  {
-    VectorBase::operator = (v);
+    template <typename number>
+    inline
+    Vector &
+    Vector::operator = (const ::Vector<number> &v)
+    {
+      VectorBase::operator = (v);
 
-    return *this;
+      return *this;
+    }
+    
   }
   
 }
index 7cc8a18451b0b5ef48fcc4741ca8925896779a6b..93371a56891cf8c2437898c2c47b3cd0558cc975 100644 (file)
 
 namespace PETScWrappers
 {
-
-
-  ParallelVector::ParallelVector ()
+  namespace MPI
   {
-                                     // this is an invalid empty vector, so we
-                                     // can just as well create a sequential
-                                     // one to avoid all the overhead incurred
-                                     // by parallelism
-    const int n = 0;
-    const int ierr
-      = VecCreateSeq (PETSC_COMM_SELF, n, &vector);
-    AssertThrow (ierr == 0, ExcPETScError(ierr));
-  }
 
-
-
-  ParallelVector::ParallelVector (const unsigned int n,
-                                  const unsigned int local_size,
-                                  const MPI_Comm    &communicator)
-                  :
-                  communicator (communicator)
-  {
-    ParallelVector::create_vector (n, local_size);
-  }
+    Vector::Vector ()
+    {
+                                       // this is an invalid empty vector, so we
+                                       // can just as well create a sequential
+                                       // one to avoid all the overhead incurred
+                                       // by parallelism
+      const int n = 0;
+      const int ierr
+        = VecCreateSeq (PETSC_COMM_SELF, n, &vector);
+      AssertThrow (ierr == 0, ExcPETScError(ierr));
+    }
+
+
+
+    Vector::Vector (const unsigned int n,
+                    const unsigned int local_size,
+                    const MPI_Comm    &communicator)
+                    :
+                    communicator (communicator)
+    {
+      Vector::create_vector (n, local_size);
+    }
 
   
 
-  ParallelVector::ParallelVector (const VectorBase &v,
-                                  const unsigned int local_size,
-                                  const MPI_Comm    &communicator)
-                  :
-                  communicator (communicator)
-  {
-    ParallelVector::create_vector (v.size(), local_size);
+    Vector::Vector (const VectorBase &v,
+                    const unsigned int local_size,
+                    const MPI_Comm    &communicator)
+                    :
+                    communicator (communicator)
+    {
+      Vector::create_vector (v.size(), local_size);
 
-    VectorBase::operator = (v);
-  }
+      VectorBase::operator = (v);
+    }
 
   
-  void
-  ParallelVector::create_vector (const unsigned int  n,
-                                 const unsigned int  local_size)
-  {
-    Assert (local_size < n, ExcIndexRange (local_size, 0, n));
+    void
+    Vector::create_vector (const unsigned int  n,
+                           const unsigned int  local_size)
+    {
+      Assert (local_size < n, ExcIndexRange (local_size, 0, n));
+
+      const int ierr
+        = VecCreateMPI (PETSC_COMM_SELF, local_size, n, &vector);
+      AssertThrow (ierr == 0, ExcPETScError(ierr));
+    }
 
-    const int ierr
-      = VecCreateMPI (PETSC_COMM_SELF, local_size, n, &vector);
-    AssertThrow (ierr == 0, ExcPETScError(ierr));
   }
+
 }
 
 #else
 // On gcc2.95 on Alpha OSF1, the native assembler does not like empty
 // files, so provide some dummy code
-namespace { void dummy () {} }
+  namespace { void dummy () {} }
 #endif // DEAL_II_USE_PETSC

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.