]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add parallel vectors.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 23 Mar 2004 19:55:56 +0000 (19:55 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 23 Mar 2004 19:55:56 +0000 (19:55 +0000)
git-svn-id: https://svn.dealii.org/trunk@8851 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/petsc_parallel_vector.h [new file with mode: 0644]
deal.II/lac/source/petsc_parallel_vector.cc [new file with mode: 0644]

diff --git a/deal.II/lac/include/lac/petsc_parallel_vector.h b/deal.II/lac/include/lac/petsc_parallel_vector.h
new file mode 100644 (file)
index 0000000..9faf31e
--- /dev/null
@@ -0,0 +1,177 @@
+//----------------------------  petsc_vector.h  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2004 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  petsc_vector.h  ---------------------------
+#ifndef __deal2__petsc_vector_h
+#define __deal2__petsc_vector_h
+
+#include <base/config.h>
+#include <base/exceptions.h>
+#include <base/subscriptor.h>
+
+#include <lac/vector.h>
+
+#ifdef DEAL_II_USE_PETSC
+
+#include <lac/petsc_vector_base.h>
+
+
+namespace PETScWrappers
+{
+/**
+ * Implementation of a parallel vector class based on PETSC and using MPI
+ * communication to synchronise distributed operations. All the functionality
+ * is actually in the base class, except for the calls to generate a
+ * parallel vector. This is possible since PETSc only works on an abstract
+ * vector type and internally distributes to functions that do the actual work
+ * depending on the actual vector type (much like using virtual
+ * functions). Only the functions creating a vector of specific type differ,
+ * and are implemented in this particular class.
+ *
+ * @author Wolfgang Bangerth, 2004
+ */
+  class ParallelVector : public VectorBase
+  {
+    public:
+                                       /**
+                                        * Default constructor. Initialize the
+                                        * vector as empty.
+                                        */
+      ParallelVector ();
+      
+                                       /**
+                                        * 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);
+    
+                                       /**
+                                        * 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);
+      
+    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);
+
+    VectorBase::operator = (v);
+  }
+
+  
+  
+  inline
+  ParallelVector &
+  ParallelVector::operator = (const PetscScalar s)
+  {
+    VectorBase::operator = (s);
+
+    return *this;
+  }
+  
+
+  template <typename number>
+  inline
+  ParallelVector &
+  ParallelVector::operator = (const ::Vector<number> &v)
+  {
+    VectorBase::operator = (v);
+
+    return *this;
+  }
+  
+}
+
+
+#endif // DEAL_II_USE_PETSC
+
+/*----------------------------   petsc_vector.h     ---------------------------*/
+
+#endif
+/*----------------------------   petsc_vector.h     ---------------------------*/
diff --git a/deal.II/lac/source/petsc_parallel_vector.cc b/deal.II/lac/source/petsc_parallel_vector.cc
new file mode 100644 (file)
index 0000000..7cc8a18
--- /dev/null
@@ -0,0 +1,79 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2004 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------------------
+
+
+#include <lac/petsc_parallel_vector.h>
+
+#include <cmath>
+
+#ifdef DEAL_II_USE_PETSC
+
+
+namespace PETScWrappers
+{
+
+
+  ParallelVector::ParallelVector ()
+  {
+                                     // 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);
+  }
+
+  
+
+  ParallelVector::ParallelVector (const VectorBase &v,
+                                  const unsigned int local_size,
+                                  const MPI_Comm    &communicator)
+                  :
+                  communicator (communicator)
+  {
+    ParallelVector::create_vector (v.size(), local_size);
+
+    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));
+
+    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 () {} }
+#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.