]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Augment VectorSpaceVector interface.
authorBruno Turcksin <bruno.turcksin@gmail.com>
Fri, 10 Feb 2017 14:26:33 +0000 (09:26 -0500)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Mon, 17 Apr 2017 12:09:52 +0000 (08:09 -0400)
14 files changed:
include/deal.II/lac/cuda_vector.h
include/deal.II/lac/la_parallel_block_vector.h
include/deal.II/lac/la_parallel_block_vector.templates.h
include/deal.II/lac/la_parallel_vector.h
include/deal.II/lac/la_parallel_vector.templates.h
include/deal.II/lac/la_vector.h
include/deal.II/lac/la_vector.templates.h
include/deal.II/lac/read_write_vector.h
include/deal.II/lac/trilinos_epetra_vector.h
include/deal.II/lac/trilinos_sparse_matrix.h
include/deal.II/lac/vector_space_vector.h
source/lac/cuda_vector.cu
source/lac/trilinos_epetra_vector.cc
tests/cuda/cuda_vector_01.cc

index 17b57dc28df8f90185e08868ea7dd04776343aba..c10ac4129ff06954710fc90d032871ccf6c99e56 100644 (file)
@@ -83,6 +83,13 @@ namespace LinearAlgebra
       void reinit(const size_type n,
                   const bool      omit_zeroing_entries = false);
 
+      /**
+       * Change the dimension to that of the vector V. The elements of V are not
+       * copied.
+       */
+      virtual void reinit(const VectorSpaceVector<Number> &V,
+                          const bool omit_zeroing_entries = false) override;
+
       /**
        * Import all the element from the input vector @p V.
        * VectorOperation::values @p operation is used to decide if the
@@ -96,6 +103,12 @@ namespace LinearAlgebra
                           std::shared_ptr<const CommunicationPatternBase> communication_pattern =
                             std::shared_ptr<const CommunicationPatternBase> ()) override;
 
+      /**
+       * Sets all elements of the vector to the scalar @p s. This operation is
+       * only allowed if @p s is equal to zero.
+       */
+      virtual Vector<Number> &operator= (const Number s) override;
+
       /**
        * Multiply the entive vector by a fixed factor.
        */
index 7bc01356919d05ae09e6b2a52c6c40da3afc71fb..7dbccdf7bd392765030f3b498fb9605db5713e97 100644 (file)
@@ -405,6 +405,13 @@ namespace LinearAlgebra
        */
       //@{
 
+      /**
+       * Change the dimension to that of the vector V. The elements of V are not
+       * copied.
+       */
+      virtual void reinit(const VectorSpaceVector<Number> &V,
+                          const bool omit_zeroing_entries = false);
+
       /**
        * Multiply the entire vector by a fixed factor.
        */
index cdce69c73a8d3daf807849caee3b990a539c2bf5..bb7e08add4380231a93e8a47cada6e287a9056db 100644 (file)
@@ -290,6 +290,18 @@ namespace LinearAlgebra
 
 
 
+    template <typename Number>
+    void
+    BlockVector<Number>::reinit(const VectorSpaceVector<Number> &V,
+                                const bool omit_zeroing_entries)
+    {
+      Assert(dynamic_cast<const BlockVector<Number> *>(&V)!=NULL,
+             ExcVectorTypeNotCompatible());
+      const BlockVector<Number> &down_V = dynamic_cast<const BlockVector<Number> &>(V);
+      reinit(down_V, omit_zeroing_entries);
+    }
+
+
     template <typename Number>
     BlockVector<Number> &
     BlockVector<Number>::operator *= (const Number factor)
index f75a71332e3a129ebe6b13cfce81c1a8bcaeb7c9..5c135c911934bc7de27e65165daf799891f26af8 100644 (file)
@@ -534,6 +534,13 @@ namespace LinearAlgebra
        */
       //@{
 
+      /**
+       * Change the dimension to that of the vector V. The elements of V are not
+       * copied.
+       */
+      virtual void reinit(const VectorSpaceVector<Number> &V,
+                          const bool omit_zeroing_entries = false);
+
       /**
        * Multiply the entire vector by a fixed factor.
        */
index fedc2724672304d2cae958095ee71ce1535ab1c2..a10d0e8a33931225a45dcd0be582ec5160ada3af 100644 (file)
@@ -939,6 +939,20 @@ namespace LinearAlgebra
 
 
 
+    template <typename Number>
+    void Vector<Number>::reinit(const VectorSpaceVector<Number> &V,
+                                const bool omit_zeroing_entries)
+    {
+      // Downcast. Throws an exception if invalid.
+      Assert(dynamic_cast<const Vector<Number> *>(&V)!=NULL,
+             ExcVectorTypeNotCompatible());
+      const Vector<Number> &down_V = dynamic_cast<const Vector<Number> &>(V);
+
+      reinit(down_V, omit_zeroing_entries);
+    }
+
+
+
     template <typename Number>
     Vector<Number> &
     Vector<Number>::operator += (const VectorSpaceVector<Number> &vv)
index fb64b714e85cfaefda296ca3926df7b8748872d2..609207b07edfff9a867295c1264aa73af909ac10 100644 (file)
@@ -87,6 +87,49 @@ namespace LinearAlgebra
     template <typename InputIterator>
     Vector(const InputIterator first, const InputIterator last);
 
+    /**
+     * Set the global size of the vector to @p size. The stored elements have
+     * their index in [0,size).
+     *
+     * If the flag @p omit_zeroing_entries is set to false, the memory will be
+     * initialized with zero, otherwise the memory will be untouched (and the
+     * user must make sure to fill it with reasonable data before using it).
+     */
+    virtual void reinit(const size_type size,
+                        const bool      omit_zeroing_entries = false);
+
+    /**
+     * Uses the same IndexSet as the one of the input vector @p in_vector and
+     * allocates memory for this vector.
+     *
+     * If the flag @p omit_zeroing_entries is set to false, the memory will be
+     * initialized with zero, otherwise the memory will be untouched (and the
+     * user must make sure to fill it with reasonable data before using it).
+     */
+    template <typename Number2>
+    void reinit(const ReadWriteVector<Number2> &in_vector,
+                const bool                      omit_zeroing_entries = false);
+
+    /**
+     * Initializes the vector. The indices are specified by @p
+     * locally_stored_indices.
+     *
+     * If the flag @p omit_zeroing_entries is set to false, the memory will be
+     * initialized with zero, otherwise the memory will be untouched (and the
+     * user must make sure to fill it with reasonable data before using it).
+     * locally_stored_indices.
+     */
+    virtual void reinit(const IndexSet &locally_stored_indices,
+                        const bool      omit_zeroing_entries = false);
+
+
+    /**
+     * Change the dimension to that of the vector V. The elements of V are not
+     * copied.
+     */
+    virtual void reinit(const VectorSpaceVector<Number> &V,
+                        const bool omit_zeroing_entries = false);
+
     /**
      * Copies the data of the input vector @p in_vector.
      */
index 84ee0877546847076a5d871b4a78a0ed2b13eeb4..c823b74777bb374d8556e4589a8c79f24965a138 100644 (file)
@@ -25,6 +25,51 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace LinearAlgebra
 {
+  template <typename Number>
+  void Vector<Number>::reinit(const size_type size,
+                              const bool omit_zeroing_entries)
+  {
+    ReadWriteVector<Number>::reinit(size, omit_zeroing_entries);
+  }
+
+
+
+  template <typename Number>
+  template <typename Number2>
+  void Vector<Number>::reinit(const ReadWriteVector<Number2> &in_vector,
+                              const bool                      omit_zeroing_entries)
+  {
+    ReadWriteVector<Number>::reinit(in_vector, omit_zeroing_entries);
+  }
+
+
+
+  template <typename Number>
+  void Vector<Number>::reinit(const IndexSet &locally_stored_indices,
+                              const bool      omit_zeroing_entries)
+  {
+    ReadWriteVector<Number>::reinit(locally_stored_indices, omit_zeroing_entries);
+  }
+
+
+
+  template <typename Number>
+  void Vector<Number>::reinit(const VectorSpaceVector<Number> &V,
+                              const bool omit_zeroing_entries)
+  {
+    // Check that casting will work.
+    Assert(dynamic_cast<const Vector<Number>*>(&V)!=NULL, ExcVectorTypeNotCompatible());
+
+    // Downcast V. If fails, throws an exception.
+    const Vector<Number> &down_V = dynamic_cast<const Vector<Number>&>(V);
+    Assert(down_V.size()==this->size(),
+           ExcMessage("Cannot add two vectors with different numbers of elements"));
+
+    ReadWriteVector<Number>::reinit(down_V, omit_zeroing_entries);
+  }
+
+
+
   template <typename Number>
   Vector<Number> &Vector<Number>::operator= (const Vector<Number> &in_vector)
   {
@@ -49,7 +94,7 @@ namespace LinearAlgebra
   {
     this->thread_loop_partitioner = in_vector.thread_loop_partitioner;
     if (this->size() != in_vector.size())
-      this->reinit(in_vector, true);
+      ReadWriteVector<Number>::reinit(in_vector, true);
 
     dealii::internal::VectorOperations::Vector_copy<Number, Number2> copier(in_vector.val, this->val);
     internal::VectorOperations::parallel_for(copier, 0, this->size(), this->thread_loop_partitioner);
@@ -446,7 +491,7 @@ namespace LinearAlgebra
 
     in.getline(buf,16,'\n');
     sz = std::atoi(buf);
-    this->reinit(sz,true);
+    ReadWriteVector<Number>::reinit(sz,true);
 
     char c;
     // in >> c;
index 2f38a41ea4da85efd38470e20d5052945c2ab853..4ae7b29fd35b5150a7f444043c218e4b1802c227 100644 (file)
@@ -31,6 +31,7 @@
 
 #ifdef DEAL_II_WITH_TRILINOS
 #include <deal.II/lac/trilinos_epetra_communication_pattern.h>
+#include <deal.II/lac/trilinos_epetra_vector.h>
 #include "Epetra_MultiVector.h"
 #endif
 
@@ -63,14 +64,6 @@ namespace TrilinosWrappers
     class Vector;
   }
 }
-
-namespace LinearAlgebra
-{
-  namespace EpetraWrappers
-  {
-    class Vector;
-  }
-}
 #endif
 
 #ifdef DEAL_II_WITH_CUDA
@@ -182,8 +175,8 @@ namespace LinearAlgebra
      * initialized with zero, otherwise the memory will be untouched (and the
      * user must make sure to fill it with reasonable data before using it).
      */
-    void reinit (const size_type size,
-                 const bool      omit_zeroing_entries = false);
+    virtual void reinit (const size_type size,
+                         const bool      omit_zeroing_entries = false);
 
     /**
      * Uses the same IndexSet as the one of the input vector @p in_vector and
@@ -206,8 +199,8 @@ namespace LinearAlgebra
      * user must make sure to fill it with reasonable data before using it).
      * locally_stored_indices.
      */
-    void reinit (const IndexSet &locally_stored_indices,
-                 const bool      omit_zeroing_entries = false);
+    virtual void reinit (const IndexSet &locally_stored_indices,
+                         const bool      omit_zeroing_entries = false);
 
 #ifdef DEAL_II_WITH_CXX11
     /**
index c08d843ca61c9becc9d2d7a623840a87d7ce07d8..1745f94768d042cb84171b27105b33538f8c8e66 100644 (file)
@@ -24,7 +24,8 @@
 #ifdef DEAL_II_WITH_MPI
 
 #include <deal.II/base/index_set.h>
-#include <deal.II/lac/read_write_vector.h>
+#include <deal.II/base/subscriptor.h>
+#include <deal.II/lac/trilinos_epetra_communication_pattern.h>
 #include <deal.II/lac/vector_space_vector.h>
 #include <memory>
 #include "mpi.h"
@@ -50,7 +51,7 @@ namespace LinearAlgebra
      * @ingroup Vectors
      * @author Bruno Turcksin, 2015
      */
-    class Vector : public VectorSpaceVector<double>
+    class Vector : public VectorSpaceVector<double>, public Subscriptor
     {
     public:
       /**
@@ -83,6 +84,13 @@ namespace LinearAlgebra
                    const MPI_Comm &communicator,
                    const bool      omit_zeroing_entries = false);
 
+      /**
+       * Change the dimension to that of the vector V. The elements of V are not
+       * copied.
+       */
+      virtual void reinit(const VectorSpaceVector<double> &V,
+                          const bool omit_zeroing_entries = false);
+
       /**
        * Copy function. This function takes a Vector and copies all the
        * elements. The Vector will have the same parallel distribution as @p
@@ -90,6 +98,12 @@ namespace LinearAlgebra
        */
       Vector &operator= (const Vector &V);
 
+      /**
+       * Sets all elements of the vector to the scalar @p s. This operation is
+       * only allowed if @p s is equal to zero.
+       */
+      Vector &operator= (const double s);
+
       /**
        * Imports all the elements present in the vector's IndexSet from the input
        * vector @p V. VectorOperation::values @p operation is used to decide if
@@ -211,6 +225,11 @@ namespace LinearAlgebra
       virtual double add_and_dot(const double a,
                                  const VectorSpaceVector<double> &V,
                                  const VectorSpaceVector<double> &W);
+      /**
+       * This function always returns false and is present only for backward
+       * compatibily.
+       */
+      bool has_ghost_elements() const;
 
       /**
        * Return the global size of the vector, equal to the sum of the number of
@@ -308,9 +327,25 @@ namespace LinearAlgebra
        */
       std::shared_ptr<const CommunicationPattern> epetra_comm_pattern;
     };
+
+
+    inline
+    bool Vector::has_ghost_elements() const
+    {
+      return false;
+    }
   }
 }
 
+
+/**
+ * Declare dealii::LinearAlgebra::EpetraWrappers::Vector as distributed vector.
+ */
+template <>
+struct is_serial_vector<LinearAlgebra::EpetraWrappers::Vector> : std_cxx11::false_type
+{
+};
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
index 15c53e6f601752c4b87ac332356187369e80e024..36cb54b780829ca3d8a5f90021c6b2745d15391b 100644 (file)
@@ -26,6 +26,7 @@
 #  include <deal.II/lac/full_matrix.h>
 #  include <deal.II/lac/exceptions.h>
 #  include <deal.II/lac/trilinos_vector.h>
+#  include <deal.II/lac/trilinos_epetra_vector.h>
 #  include <deal.II/lac/vector_view.h>
 
 #ifdef DEAL_II_WITH_CXX11
index ee966cc76da5ef90caf0d2c9f1be3d4f5a2a4445..8a18ff9f657e379c219f6bc03b10634030939b9d 100644 (file)
@@ -55,6 +55,18 @@ namespace LinearAlgebra
     typedef types::global_dof_index                           size_type;
     typedef typename numbers::NumberTraits<Number>::real_type real_type;
 
+    /**
+     * Change the dimension to that of the vector V. The elements of V are not
+     * copied.
+     */
+    virtual void reinit(const VectorSpaceVector<Number> &V,
+                        const bool omit_zeroing_entries = false) = 0;
+
+    /**
+     * Sets all elements of the vector to the scalar @p s. This operation is
+     * only allowed if @p s is equal to zero.
+     */
+    virtual VectorSpaceVector<Number> &operator= (const Number s) = 0;
 
     /**
      * Multiply the entire vector by a fixed factor.
index 64e60a6a8c5b814804bb35cf7aecba16ec3313e1..4d0921f2f0d6f6c0ec35c4cf081c5b5b1e8eca5c 100644 (file)
@@ -564,6 +564,15 @@ namespace LinearAlgebra
 
 
 
+    template <typename Number>
+    void Vector<Number>::reinit(const VectorSpaceVector<Number> &V,
+                                const bool omit_zeroing_entries)
+    {
+      reinit(V.size(), omit_zeroing_entries);
+    }
+
+
+
     template <typename Number>
     void Vector<Number>::import(const ReadWriteVector<Number> &V,
                                 VectorOperation::values operation,
@@ -606,6 +615,20 @@ namespace LinearAlgebra
 
 
 
+    template <typename Number>
+    Vector<Number> &Vector<Number>::operator= (const Number s)
+    {
+      Assert(s == Number(), ExcMessage("Onlyt 0 can be assigned to a vector."));
+      (void)s;
+
+      cudaError_t error_code = cudaMemset(val, 0, n_elements*sizeof(Number));
+      AssertCuda(error_code);
+
+      return *this;
+    }
+
+
+
     template <typename Number>
     Vector<Number> &Vector<Number>::operator*= (const Number factor)
     {
@@ -1029,7 +1052,7 @@ namespace LinearAlgebra
     void Vector<Number>::print(std::ostream       &out,
                                const unsigned int  precision,
                                const bool          scientific,
-                               const bool          across) const
+                               const bool          ) const
     {
       AssertThrow(out, ExcIO());
       std::ios::fmtflags old_flags = out.flags();
index 078c36acf320f89e8eb9823124ccbd18540c3216..e60297315a0ccdb20fe0bf61d837f994c5d1f792 100644 (file)
@@ -23,6 +23,7 @@
 
 #include <boost/io/ios_state.hpp>
 
+#include <deal.II/lac/read_write_vector.h>
 #include "Epetra_Import.h"
 #include "Epetra_Map.h"
 #include "Epetra_MpiComm.h"
@@ -43,6 +44,7 @@ namespace LinearAlgebra
 
     Vector::Vector(const Vector &V)
       :
+      Subscriptor(),
       vector(new Epetra_FEVector(V.trilinos_vector()))
     {}
 
@@ -73,6 +75,21 @@ namespace LinearAlgebra
 
 
 
+    void Vector::reinit(const VectorSpaceVector<double> &V,
+                        const bool omit_zeroing_entries)
+    {
+      // Check that casting will work.
+      Assert(dynamic_cast<const Vector *>(&V)!=NULL, ExcVectorTypeNotCompatible());
+
+      // Downcast V. If fails, throws an exception.
+      const Vector &down_V = dynamic_cast<const Vector &>(V);
+
+      reinit(down_V.locally_owned_elements(), down_V.get_mpi_communicator(),
+             omit_zeroing_entries);
+    }
+
+
+
     Vector &Vector::operator= (const Vector &V)
     {
       // Distinguish three cases:
@@ -100,6 +117,16 @@ namespace LinearAlgebra
 
 
 
+    Vector &Vector::operator= (const double s)
+    {
+      Assert(s==0., ExcMessage("Only 0 can be assigned to a vector."));
+      vector->PutScalar(s);
+
+      return *this;
+    }
+
+
+
     void Vector::import(const ReadWriteVector<double>                  &V,
                         VectorOperation::values                         operation,
                         std::shared_ptr<const CommunicationPatternBase> communication_pattern)
index f3c8b53685a6435ccf226c997db8751e24ac958d..2632adc84aa3ad90094301cf74a769566c550504 100644 (file)
@@ -29,10 +29,13 @@ void test()
   LinearAlgebra::CUDAWrappers::Vector<double> a;
   LinearAlgebra::CUDAWrappers::Vector<double> b(size);
   LinearAlgebra::CUDAWrappers::Vector<double> c(b);
+  LinearAlgebra::CUDAWrappers::Vector<double> d;
+  d.reinit(c);
 
   AssertThrow(a.size()==0, ExcMessage("Vector has the wrong size."));
   AssertThrow(b.size()==size, ExcMessage("Vector has the wrong size."));
   AssertThrow(c.size()==size, ExcMessage("Vector has the wrong size."));
+  AssertThrow(d.size()==size, ExcMessage("Vector has the wrong size."));
 
   a.reinit(size);
   AssertThrow(a.size()==size, ExcMessage("Vector has the wrong size."));
@@ -95,6 +98,11 @@ void test()
   c.import(read_write_1, VectorOperation::insert);
   const double val = b*c;
   AssertThrow(val==328350., ExcMessage("Problem in operator *."));
+
+  b = 0.;
+  read_write_3.import(b, VectorOperation::insert);
+  for (unsigned int i=0; i<size; ++i)
+    AssertThrow(read_write_3[i] == 0.,ExcMessage("Problem in operator =."));
 }
 
 int main(int argc, char **argv)

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.