]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Make parallel vector derive from VectorSpaceVector
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 18 May 2016 15:46:32 +0000 (17:46 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 14 Jun 2016 15:10:11 +0000 (17:10 +0200)
20 files changed:
include/deal.II/base/partitioner.h
include/deal.II/lac/la_parallel_vector.h [new file with mode: 0644]
include/deal.II/lac/la_parallel_vector.templates.h [moved from include/deal.II/lac/parallel_vector.templates.h with 50% similarity]
include/deal.II/lac/parallel_vector.h
include/deal.II/lac/precondition.h
include/deal.II/lac/read_write_vector.h
include/deal.II/lac/vector_operations_internal.h
include/deal.II/lac/vector_space_vector.h
include/deal.II/matrix_free/fe_evaluation.h
source/base/partitioner.cc
source/lac/CMakeLists.txt
source/lac/la_parallel_vector.cc [new file with mode: 0644]
source/lac/la_parallel_vector.inst.in [moved from source/lac/parallel_vector.inst.in with 90% similarity]
source/lac/parallel_vector.cc [deleted file]
tests/mpi/parallel_vector_06.cc
tests/mpi/parallel_vector_06.mpirun=1.output
tests/mpi/parallel_vector_06.mpirun=10.output
tests/mpi/parallel_vector_06.mpirun=4.output
tests/mpi/parallel_vector_08.cc
tests/mpi/parallel_vector_11.cc

index 415782de2def42cb2cb35ea9986ee50edf07350c..0ae382e1782425c54ee3fa2c3c2e2f340476deca 100644 (file)
@@ -22,6 +22,7 @@
 #include <deal.II/base/types.h>
 #include <deal.II/base/utilities.h>
 #include <deal.II/base/memory_consumption.h>
+#include <deal.II/lac/communication_pattern_base.h>
 
 #include <limits>
 
@@ -58,7 +59,7 @@ namespace Utilities
      *
      * @author Katharina Kormann, Martin Kronbichler, 2010, 2011
      */
-    class Partitioner
+    class Partitioner : public ::dealii::LinearAlgebra::CommunicationPatternBase
     {
     public:
       /**
@@ -93,6 +94,17 @@ namespace Utilities
       Partitioner (const IndexSet &locally_owned_indices,
                    const MPI_Comm  communicator_in);
 
+      /**
+       * Reinitialize the communication pattern. The first argument @p
+       * vector_space_vector_index_set is the index set associated to a
+       * VectorSpaceVector object. The second argument @p
+       * read_write_vector_index_set is the index set associated to a
+       * ReadWriteVector object.
+       */
+      virtual void reinit(const IndexSet &vector_space_vector_index_set,
+                          const IndexSet &read_write_vector_index_set,
+                          const MPI_Comm &communicator);
+
       /**
        * Sets the locally owned indices. Used in the constructor.
        */
@@ -250,7 +262,12 @@ namespace Utilities
       /**
        * Returns the MPI communicator underlying the partitioner object.
        */
-      const MPI_Comm &get_communicator() const;
+      const MPI_Comm &get_communicator() const DEAL_II_DEPRECATED;
+
+      /**
+       * Returns the MPI communicator underlying the partitioner object.
+       */
+      virtual const MPI_Comm &get_mpi_communicator() const;
 
       /**
        * Returns whether ghost indices have been explicitly added as a @p
@@ -277,7 +294,7 @@ namespace Utilities
       /**
        * The global size of the vector over all processors
        */
-      const types::global_dof_index global_size;
+      types::global_dof_index global_size;
 
       /**
        * The range of the vector that is stored locally.
@@ -341,7 +358,7 @@ namespace Utilities
       /**
        * The MPI communicator involved in the problem
        */
-      const MPI_Comm communicator;
+      MPI_Comm communicator;
 
       /**
        * Stores whether the ghost indices have been explicitly set.
@@ -533,6 +550,15 @@ namespace Utilities
 
 
 
+    inline
+    const MPI_Comm &
+    Partitioner::get_mpi_communicator() const
+    {
+      return communicator;
+    }
+
+
+
     inline
     bool
     Partitioner::ghost_indices_initialized() const
diff --git a/include/deal.II/lac/la_parallel_vector.h b/include/deal.II/lac/la_parallel_vector.h
new file mode 100644 (file)
index 0000000..e9bf996
--- /dev/null
@@ -0,0 +1,1383 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2011 - 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii__la_parallel_vector_h
+#define dealii__la_parallel_vector_h
+
+#include <deal.II/base/config.h>
+#include <deal.II/base/mpi.h>
+#include <deal.II/base/numbers.h>
+#include <deal.II/base/partitioner.h>
+#include <deal.II/base/std_cxx11/shared_ptr.h>
+#include <deal.II/base/thread_management.h>
+#include <deal.II/lac/vector_space_vector.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace parallel
+{
+  namespace distributed
+  {
+    template <typename> class BlockVector;
+  }
+}
+
+namespace LinearAlgebra
+{
+  template <typename> class ReadWriteVector;
+}
+
+namespace LinearAlgebra
+{
+  namespace parallel
+  {
+    /*! @addtogroup Vectors
+     *@{
+     */
+
+    /**
+     * Implementation of a parallel vector class. The design of this class is
+     * similar to the standard ::dealii::Vector class in deal.II, with the
+     * exception that storage is distributed with MPI.
+     *
+     * The vector is designed for the following scheme of parallel
+     * partitioning:
+     * <ul>
+     * <li> The indices held by individual processes (locally owned part) in
+     * the MPI parallelization form a contiguous range
+     * <code>[my_first_index,my_last_index)</code>.
+     * <li> Ghost indices residing on arbitrary positions of other processors
+     * are allowed. It is in general more efficient if ghost indices are
+     * clustered, since they are stored as a set of intervals. The
+     * communication pattern of the ghost indices is determined when calling
+     * the function <code>reinit (locally_owned, ghost_indices,
+     * communicator)</code>, and retained until the partitioning is changed.
+     * This allows for efficient parallel communication of indices. In
+     * particular, it stores the communication pattern, rather than having to
+     * compute it again for every communication. For more information on ghost
+     * vectors, see also the
+     * @ref GlossGhostedVector "glossary entry on vectors with ghost elements".
+     * <li> Besides the usual global access operator() it is also possible to
+     * access vector entries in the local index space with the function @p
+     * local_element(). Locally owned indices are placed first, [0,
+     * local_size()), and then all ghost indices follow after them
+     * contiguously, [local_size(), local_size()+n_ghost_entries()).
+     * </ul>
+     *
+     * Functions related to parallel functionality:
+     * <ul>
+     * <li> The function <code>compress()</code> goes through the data
+     * associated with ghost indices and communicates it to the owner process,
+     * which can then add it to the correct position. This can be used e.g.
+     * after having run an assembly routine involving ghosts that fill this
+     * vector. Note that the @p insert mode of @p compress() does not set the
+     * elements included in ghost entries but simply discards them, assuming
+     * that the owning processor has set them to the desired value already
+     * (See also the
+     * @ref GlossCompress "glossary entry on compress").
+     * <li> The <code>update_ghost_values()</code> function imports the data
+     * from the owning processor to the ghost indices in order to provide read
+     * access to the data associated with ghosts.
+     * <li> It is possible to split the above functions into two phases, where
+     * the first initiates the communication and the second one finishes it.
+     * These functions can be used to overlap communication with computations
+     * in other parts of the code.
+     * <li> Of course, reduction operations (like norms) make use of
+     * collective all-to-all MPI communications.
+     * </ul>
+     *
+     * This vector can take two different states with respect to ghost
+     * elements:
+     * <ul>
+     * <li> After creation and whenever zero_out_ghosts() is called (or
+     * <code>operator= (0.)</code>), the vector does only allow writing into
+     * ghost elements but not reading from ghost elements.
+     * <li> After a call to update_ghost_values(), the vector does not allow
+     * writing into ghost elements but only reading from them. This is to
+     * avoid undesired ghost data artifacts when calling compress() after
+     * modifying some vector entries. The current status of the ghost entries
+     * (read mode or write mode) can be queried by the method
+     * has_ghost_elements(), which returns <code>true</code> exactly when
+     * ghost elements have been updated and <code>false</code> otherwise,
+     * irrespective of the actual number of ghost entries in the vector layout
+     * (for that information, use n_ghost_entries() instead).
+     * </ul>
+     *
+     * This vector uses the facilities of the class dealii::Vector<Number> for
+     * implementing the operations on the local range of the vector. In
+     * particular, it also inherits thread parallelism that splits most
+     * vector-vector operations into smaller chunks if the program uses
+     * multiple threads. This may or may not be desired when working also with
+     * MPI.
+     *
+     * <h4>Limitations regarding the vector size</h4>
+     *
+     * This vector class is based on two different number types for indexing.
+     * The so-called global index type encodes the overall size of the vector.
+     * Its type is types::global_dof_index. The largest possible value is
+     * <code>2^32-1</code> or approximately 4 billion in case 64 bit integers
+     * are disabled at configuration of deal.II (default case) or
+     * <code>2^64-1</code> or approximately <code>10^19</code> if 64 bit
+     * integers are enabled (see the glossary entry on
+     * @ref GlobalDoFIndex
+     * for further information).
+     *
+     * The second relevant index type is the local index used within one MPI
+     * rank. As opposed to the global index, the implementation assumes 32-bit
+     * unsigned integers unconditionally. In other words, to actually use a
+     * vector with more than four billion entries, you need to use MPI with
+     * more than one rank (which in general is a safe assumption since four
+     * billion entries consume at least 16 GB of memory for floats or 32 GB of
+     * memory for doubles) and enable 64-bit indices. If more than 4 billion
+     * local elements are present, the implementation tries to detect that,
+     * which triggers an exception and aborts the code. Note, however, that
+     * the detection of overflow is tricky and the detection mechanism might
+     * fail in some circumstances. Therefore, it is strongly recommended to
+     * not rely on this class to automatically detect the unsupported case.
+     *
+     * @author Katharina Kormann, Martin Kronbichler, 2010, 2011, 2016
+     */
+    template <typename Number>
+    class Vector : public ::dealii::LinearAlgebra::VectorSpaceVector<Number>,
+      public Subscriptor
+    {
+    public:
+      typedef Number                                            value_type;
+      typedef value_type                                       *pointer;
+      typedef const value_type                                 *const_pointer;
+      typedef value_type                                       *iterator;
+      typedef const value_type                                 *const_iterator;
+      typedef value_type                                       &reference;
+      typedef const value_type                                 &const_reference;
+      typedef types::global_dof_index                           size_type;
+      typedef typename numbers::NumberTraits<Number>::real_type real_type;
+
+      /**
+       * A variable that indicates whether this vector supports distributed
+       * data storage. If true, then this vector also needs an appropriate
+       * compress() function that allows communicating recent set or add
+       * operations to individual elements to be communicated to other
+       * processors.
+       *
+       * For the current class, the variable equals true, since it does
+       * support parallel data storage.
+       */
+      static const bool supports_distributed_data = true;
+
+      /**
+       * @name 1: Basic Object-handling
+       */
+      //@{
+      /**
+       * Empty constructor.
+       */
+      Vector ();
+
+      /**
+       * Copy constructor. Uses the parallel partitioning of @p in_vector.
+       */
+      Vector (const Vector<Number> &in_vector);
+
+      /**
+       * Constructs a parallel vector of the given global size without any
+       * actual parallel distribution.
+       */
+      Vector (const size_type size);
+
+      /**
+       * Constructs a parallel vector. The local range is specified by @p
+       * locally_owned_set (note that this must be a contiguous interval,
+       * multiple intervals are not possible). The IndexSet @p ghost_indices
+       * specifies ghost indices, i.e., indices which one might need to read
+       * data from or accumulate data from. It is allowed that the set of
+       * ghost indices also contains the local range, but it does not need to.
+       *
+       * This function involves global communication, so it should only be
+       * called once for a given layout. Use the constructor with
+       * Vector<Number> argument to create additional vectors with the same
+       * parallel layout.
+       *
+       * @see
+       * @ref GlossGhostedVector "vectors with ghost elements"
+       */
+      Vector (const IndexSet &local_range,
+              const IndexSet &ghost_indices,
+              const MPI_Comm  communicator);
+
+      /**
+       * Same constructor as above but without any ghost indices.
+       */
+      Vector (const IndexSet &local_range,
+              const MPI_Comm  communicator);
+
+      /**
+       * Create the vector based on the parallel partitioning described in @p
+       * partitioner. The input argument is a shared pointer, which store the
+       * partitioner data only once and share it between several vectors with
+       * the same layout.
+       */
+      Vector (const std_cxx11::shared_ptr<const Utilities::MPI::Partitioner> &partitioner);
+
+      /**
+       * Destructor.
+       */
+      virtual ~Vector ();
+
+      /**
+       * Sets the global size of the vector to @p size without any actual
+       * parallel distribution.
+       */
+      void reinit (const size_type size,
+                   const bool      omit_zeroing_entries = false);
+
+      /**
+       * Uses the parallel layout of the input vector @p in_vector and
+       * allocates memory for this vector. Recommended initialization function
+       * when several vectors with the same layout should be created.
+       *
+       * 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 Vector<Number2> &in_vector,
+                  const bool             omit_zeroing_entries = false);
+
+      /**
+       * Initialize the vector. The local range is specified by @p
+       * locally_owned_set (note that this must be a contiguous interval,
+       * multiple intervals are not possible). The IndexSet @p ghost_indices
+       * specifies ghost indices, i.e., indices which one might need to read
+       * data from or accumulate data from. It is allowed that the set of
+       * ghost indices also contains the local range, but it does not need to.
+       *
+       * This function involves global communication, so it should only be
+       * called once for a given layout. Use the @p reinit function with
+       * Vector<Number> argument to create additional vectors with the same
+       * parallel layout.
+       *
+       * @see
+       * @ref GlossGhostedVector "vectors with ghost elements"
+       */
+      void reinit (const IndexSet &local_range,
+                   const IndexSet &ghost_indices,
+                   const MPI_Comm  communicator);
+
+      /**
+       * Same as above, but without ghost entries.
+       */
+      void reinit (const IndexSet &local_range,
+                   const MPI_Comm  communicator);
+
+      /**
+       * Initialize the vector given to the parallel partitioning described in
+       * @p partitioner. The input argument is a shared pointer, which store
+       * the partitioner data only once and share it between several vectors
+       * with the same layout.
+       */
+      void reinit (const std_cxx11::shared_ptr<const Utilities::MPI::Partitioner> &partitioner);
+
+      /**
+       * Swap the contents of this vector and the other vector @p v. One could
+       * do this operation with a temporary variable and copying over the data
+       * elements, but this function is significantly more efficient since it
+       * only swaps the pointers to the data of the two vectors and therefore
+       * does not need to allocate temporary storage and move data around.
+       *
+       * This function is analog to the the @p swap function of all C++
+       * standard containers. Also, there is a global function
+       * <tt>swap(u,v)</tt> that simply calls <tt>u.swap(v)</tt>, again in
+       * analogy to standard functions.
+       *
+       * This function is virtual in order to allow for derived classes to
+       * handle memory separately.
+       */
+      void swap (Vector<Number> &v);
+
+      /**
+       * Assigns the vector to the parallel partitioning of the input vector
+       * @p in_vector, and copies all the data.
+       *
+       * If one of the input vector or the calling vector (to the left of the
+       * assignment operator) had ghost elements set before this operation,
+       * the calling vector will have ghost values set. Otherwise, it will be
+       * in write mode. If the input vector does not have any ghost elements
+       * at all, the vector will also update its ghost values in analogy to
+       * the respective setting the Trilinos and PETSc vectors.
+       */
+      Vector<Number> &
+      operator = (const Vector<Number> &in_vector);
+
+      /**
+       * Assigns the vector to the parallel partitioning of the input vector
+       * @p in_vector, and copies all the data.
+       *
+       * If one of the input vector or the calling vector (to the left of the
+       * assignment operator) had ghost elements set before this operation,
+       * the calling vector will have ghost values set. Otherwise, it will be
+       * in write mode. If the input vector does not have any ghost elements
+       * at all, the vector will also update its ghost values in analogy to
+       * the respective setting the Trilinos and PETSc vectors.
+       */
+      template <typename Number2>
+      Vector<Number> &
+      operator = (const Vector<Number2> &in_vector);
+      //@}
+
+      /**
+       * @name 2: Parallel data exchange
+       */
+      //@{
+      /**
+       * This function copies the data that has accumulated in the data buffer
+       * for ghost indices to the owning processor. For the meaning of the
+       * argument @p operation, see the entry on
+       * @ref GlossCompress "Compressing distributed vectors and matrices"
+       * in the glossary.
+       *
+       * There are two variants for this function. If called with argument @p
+       * VectorOperation::add adds all the data accumulated in ghost elements
+       * to the respective elements on the owning processor and clears the
+       * ghost array afterwards. If called with argument @p
+       * VectorOperation::insert, a set operation is performed. Since setting
+       * elements in a vector with ghost elements is ambiguous (as one can set
+       * both the element on the ghost site as well as the owning site), this
+       * operation makes the assumption that all data is set correctly on the
+       * owning processor. Upon call of compress(VectorOperation::insert), all
+       * ghost entries are thus simply zeroed out (using zero_ghost_values()).
+       * In debug mode, a check is performed for whether the data set is
+       * actually consistent between processors, i.e., whenever a non-zero
+       * ghost element is found, it is compared to the value on the owning
+       * processor and an exception is thrown if these elements do not agree.
+       */
+      void compress (::dealii::VectorOperation::values operation);
+
+      /**
+       * Fills the data field for ghost indices with the values stored in the
+       * respective positions of the owning processor. This function is needed
+       * before reading from ghosts. The function is @p const even though
+       * ghost data is changed. This is needed to allow functions with a @p
+       * const vector to perform the data exchange without creating
+       * temporaries.
+       *
+       * After calling this method, write access to ghost elements of the
+       * vector is forbidden and an exception is thrown. Only read access to
+       * ghost elements is allowed in this state. Note that all subsequent
+       * operations on this vector, like global vector addition, etc., will
+       * also update the ghost values by a call to this method after the
+       * operation. However, global reduction operations like norms or the
+       * inner product will always ignore ghost elements in order to avoid
+       * counting the ghost data more than once. To allow writing to ghost
+       * elements again, call zero_out_ghosts().
+       *
+       * @see
+       * @ref GlossGhostedVector "vectors with ghost elements"
+       */
+      void update_ghost_values () const;
+
+      /**
+       * Initiates communication for the @p compress() function with non-
+       * blocking communication. This function does not wait for the transfer
+       * to finish, in order to allow for other computations during the time
+       * it takes until all data arrives.
+       *
+       * Before the data is actually exchanged, the function must be followed
+       * by a call to @p compress_finish().
+       *
+       * In case this function is called for more than one vector before @p
+       * compress_finish() is invoked, it is mandatory to specify a unique
+       * communication channel to each such call, in order to avoid several
+       * messages with the same ID that will corrupt this operation.
+       */
+      void compress_start (const unsigned int communication_channel = 0,
+                           ::dealii::VectorOperation::values operation = VectorOperation::add);
+
+      /**
+       * For all requests that have been initiated in compress_start, wait for
+       * the communication to finish. Once it is finished, add or set the data
+       * (depending on the flag operation) to the respective positions in the
+       * owning processor, and clear the contents in the ghost data fields.
+       * The meaning of this argument is the same as in compress().
+       *
+       * This function should be called exactly once per vector after calling
+       * compress_start, otherwise the result is undefined. In particular, it
+       * is not well-defined to call compress_start on the same vector again
+       * before compress_finished has been called. However, there is no
+       * warning to prevent this situation.
+       *
+       * Must follow a call to the @p compress_start function.
+       */
+      void compress_finish (::dealii::VectorOperation::values operation);
+
+      /**
+       * Initiates communication for the @p update_ghost_values() function
+       * with non-blocking communication. This function does not wait for the
+       * transfer to finish, in order to allow for other computations during
+       * the time it takes until all data arrives.
+       *
+       * Before the data is actually exchanged, the function must be followed
+       * by a call to @p update_ghost_values_finish().
+       *
+       * In case this function is called for more than one vector before @p
+       * update_ghost_values_finish() is invoked, it is mandatory to specify a
+       * unique communication channel to each such call, in order to avoid
+       * several messages with the same ID that will corrupt this operation.
+       */
+      void update_ghost_values_start (const unsigned int communication_channel = 0) const;
+
+
+      /**
+       * For all requests that have been started in update_ghost_values_start,
+       * wait for the communication to finish.
+       *
+       * Must follow a call to the @p update_ghost_values_start function
+       * before reading data from ghost indices.
+       */
+      void update_ghost_values_finish () const;
+
+      /**
+       * This method zeros the entries on ghost dofs, but does not touch
+       * locally owned DoFs.
+       *
+       * After calling this method, read access to ghost elements of the
+       * vector is forbidden and an exception is thrown. Only write access to
+       * ghost elements is allowed in this state.
+       */
+      void zero_out_ghosts ();
+
+      /**
+       * Returns whether the vector currently is in a state where ghost values
+       * can be read or not. This is the same functionality as other parallel
+       * vectors have. If this method returns false, this only means that
+       * read-access to ghost elements is prohibited whereas write access is
+       * still possible (to those entries specified as ghosts during
+       * initialization), not that there are no ghost elements at all.
+       *
+       * @see
+       * @ref GlossGhostedVector "vectors with ghost elements"
+       */
+      bool has_ghost_elements() const;
+      //@}
+
+      /**
+       * @name 3: Implementation of VectorSpaceVector
+       */
+      //@{
+
+      /**
+       * Multiply the entire vector by a fixed factor.
+       */
+      virtual Vector<Number> &operator*= (const Number factor);
+
+      /**
+       * Divide the entire vector by a fixed factor.
+       */
+      virtual Vector<Number> &operator/= (const Number factor);
+
+      /**
+       * Add the vector @p V to the present one.
+       */
+      virtual Vector<Number> &operator+= (const VectorSpaceVector<Number> &V);
+
+      /**
+       * Subtract the vector @p V from the present one.
+       */
+      virtual Vector<Number> &operator-= (const VectorSpaceVector<Number> &V);
+
+      /**
+       * Import all the elements present in the vector's IndexSet from the input
+       * vector @p V. VectorOperation::values @p operation is used to decide if
+       * the elements in @p V should be added to the current vector or replace the
+       * current elements. The last parameter can be used if the same
+       * communication pattern is used multiple times. This can be used to improve
+       * performance.
+       */
+      virtual void import(const LinearAlgebra::ReadWriteVector<Number> &V,
+                          VectorOperation::values operation,
+                          std_cxx11::shared_ptr<const CommunicationPatternBase> communication_pattern =
+                            std_cxx11::shared_ptr<const CommunicationPatternBase> ());
+
+      /**
+       * Return the scalar product of two vectors.
+       */
+      virtual Number operator* (const VectorSpaceVector<Number> &V) const;
+
+      /**
+       * Add @p a to all components. Note that @p a is a scalar not a vector.
+       */
+      virtual void add(const Number a);
+
+      /**
+       * Simple addition of a multiple of a vector, i.e. <tt>*this += a*V</tt>.
+       */
+      virtual void add(const Number a, const VectorSpaceVector<Number> &V);
+
+      /**
+       * Multiple addition of scaled vectors, i.e. <tt>*this += a*V+b*W</tt>.
+       */
+      virtual void add(const Number a, const VectorSpaceVector<Number> &V,
+                       const Number b, const VectorSpaceVector<Number> &W);
+
+      /**
+       * Scaling and simple addition of a multiple of a vector, i.e. <tt>*this =
+       * s*(*this)+a*V</tt>.
+       */
+      virtual void sadd(const Number s, const Number a,
+                        const VectorSpaceVector<Number> &V);
+
+      /**
+       * Scale each element of this vector by the corresponding element in the
+       * argument. This function is mostly meant to simulate multiplication (and
+       * immediate re-assignment) by a diagonal scaling matrix.
+       */
+      virtual void scale(const VectorSpaceVector<Number> &scaling_factors);
+
+      /**
+       * Assignment <tt>*this = a*V</tt>.
+       */
+      virtual void equ(const Number a, const VectorSpaceVector<Number> &V);
+
+      /**
+       * Return the l<sub>1</sub> norm of the vector (i.e., the sum of the
+       * absolute values of all entries among all processors).
+       */
+      virtual real_type l1_norm() const;
+
+      /**
+       * Return the l<sub>2</sub> norm of the vector (i.e., the square root of
+       * the sum of the square of all entries among all processors).
+       */
+      virtual real_type l2_norm() const;
+
+      /**
+       * Return the maximum norm of the vector (i.e., the maximum absolute value
+       * among all entries and among all processors).
+       */
+      virtual real_type linfty_norm() const;
+
+      /**
+       * Perform a combined operation of a vector addition and a subsequent
+       * inner product, returning the value of the inner product. In other
+       * words, the result of this function is the same as if the user called
+       * @code
+       * this->add(a, V);
+       * return_value = *this * W;
+       * @endcode
+       *
+       * The reason this function exists is that this operation involves less
+       * memory transfer than calling the two functions separately. This method
+       * only needs to load three vectors, @p this, @p V, @p W, whereas calling
+       * separate methods means to load the calling vector @p this twice. Since
+       * most vector operations are memory transfer limited, this reduces the
+       * time by 25\% (or 50\% if @p W equals @p this).
+       */
+      virtual Number add_and_dot(const Number a,
+                                 const VectorSpaceVector<Number> &V,
+                                 const VectorSpaceVector<Number> &W);
+
+      /**
+       * Return the global size of the vector, equal to the sum of the number of
+       * locally owned indices among all processors.
+       */
+      virtual size_type size() const;
+
+      /**
+       * Return an index set that describes which elements of this vector are
+       * owned by the current processor. As a consequence, the index sets
+       * returned on different procesors if this is a distributed vector will
+       * form disjoint sets that add up to the complete index set. Obviously, if
+       * a vector is created on only one processor, then the result would
+       * satisfy
+       * @code
+       *  vec.locally_owned_elements() == complete_index_set(vec.size())
+       * @endcode
+       */
+      virtual dealii::IndexSet locally_owned_elements() const;
+
+      /**
+       * Print the vector to the output stream @p out.
+       */
+      virtual void print(std::ostream &out,
+                         const unsigned int precision=3,
+                         const bool scientific=true,
+                         const bool across=true) const;
+
+      /**
+       * Return the memory consumption of this class in bytes.
+       */
+      virtual std::size_t memory_consumption() const;
+      //@}
+
+      /**
+       * @name 4: Other vector operations not included in VectorSpaceVector
+       */
+      //@{
+
+      /**
+       * Sets all elements of the vector to the scalar @p s. If the scalar is
+       * zero, also ghost elements are set to zero, otherwise they remain
+       * unchanged.
+       */
+      Vector<Number> &operator = (const Number s);
+
+      /**
+       * This is a collective add operation that adds a whole set of values
+       * stored in @p values to the vector components specified by @p indices.
+       */
+      template <typename OtherNumber>
+      void add (const std::vector<size_type>        &indices,
+                const ::dealii::Vector<OtherNumber> &values);
+
+      /**
+       * Scaling and simple vector addition, i.e.  <tt>*this =
+       * s*(*this)+V</tt>.
+       */
+      void sadd (const Number          s,
+                 const Vector<Number> &V);
+
+      /**
+       * Scaling and multiple addition.
+       *
+       * This function is deprecated.
+       */
+      void sadd (const Number          s,
+                 const Number          a,
+                 const Vector<Number> &V,
+                 const Number          b,
+                 const Vector<Number> &W) DEAL_II_DEPRECATED;
+
+      /**
+       * Assignment <tt>*this = a*u + b*v</tt>.
+       *
+       * This function is deprecated.
+       */
+      void equ (const Number a, const Vector<Number> &u,
+                const Number b, const Vector<Number> &v) DEAL_II_DEPRECATED;
+
+      //@}
+
+
+      /**
+       * @name 5: Entry access and local data representation
+       */
+      //@{
+
+      /**
+       * Returns the local size of the vector, i.e., the number of indices
+       * owned locally.
+       */
+      size_type local_size() const;
+
+      /**
+       * Returns the half-open interval that specifies the locally owned range
+       * of the vector. Note that <code>local_size() == local_range().second -
+       * local_range().first</code>.
+       *
+       * This function is deprecated.
+       */
+      std::pair<size_type, size_type> local_range () const DEAL_II_DEPRECATED;
+
+      /**
+       * Returns true if the given global index is in the local range of this
+       * processor.
+       *
+       * This function is deprecated.
+       */
+      bool in_local_range (const size_type global_index) const DEAL_II_DEPRECATED;
+
+      /**
+       * Returns the number of ghost elements present on the vector.
+       *
+       * This function is deprecated.
+       */
+      size_type n_ghost_entries () const DEAL_II_DEPRECATED;
+
+      /**
+       * Return an index set that describes which elements of this vector are
+       * not owned by the current processor but can be written into or read
+       * from locally (ghost elements).
+       *
+       * This function is deprecated.
+       */
+      const IndexSet &ghost_elements() const DEAL_II_DEPRECATED;
+
+      /**
+       * Returns whether the given global index is a ghost index on the
+       * present processor. Returns false for indices that are owned locally
+       * and for indices not present at all.
+       *
+       * This function is deprecated.
+       */
+      bool is_ghost_entry (const types::global_dof_index global_index) const DEAL_II_DEPRECATED;
+
+      /**
+       * Make the @p Vector class a bit like the <tt>vector<></tt> class of
+       * the C++ standard library by returning iterators to the start and end
+       * of the <i>locally owned</i> elements of this vector.
+       *
+       * It holds that end() - begin() == local_size().
+       */
+      iterator begin ();
+
+      /**
+       * Return constant iterator to the start of the locally owned elements
+       * of the vector.
+       */
+      const_iterator begin () const;
+
+      /**
+       * Return an iterator pointing to the element past the end of the array
+       * of locally owned entries.
+       */
+      iterator end ();
+
+      /**
+       * Return a constant iterator pointing to the element past the end of
+       * the array of the locally owned entries.
+       */
+      const_iterator end () const;
+
+      /**
+       * Read access to the data in the position corresponding to @p
+       * global_index. The index must be either in the local range of the
+       * vector or be specified as a ghost index at construction.
+       *
+       * Performance: <tt>O(1)</tt> for locally owned elements that represent
+       * a contiguous range and <tt>O(log(n<sub>ranges</sub>))</tt> for ghost
+       * elements (quite fast, but slower than local_element()).
+       */
+      Number operator () (const size_type global_index) const;
+
+      /**
+       * Read and write access to the data in the position corresponding to @p
+       * global_index. The index must be either in the local range of the
+       * vector or be specified as a ghost index at construction.
+       *
+       * Performance: <tt>O(1)</tt> for locally owned elements that represent
+       * a contiguous range and <tt>O(log(n<sub>ranges</sub>))</tt> for ghost
+       * elements (quite fast, but slower than local_element()).
+       */
+      Number &operator () (const size_type global_index);
+
+      /**
+       * Read access to the data in the position corresponding to @p
+       * global_index. The index must be either in the local range of the
+       * vector or be specified as a ghost index at construction.
+       *
+       * This function does the same thing as operator().
+       */
+      Number operator [] (const size_type global_index) const;
+      /**
+       * Read and write access to the data in the position corresponding to @p
+       * global_index. The index must be either in the local range of the
+       * vector or be specified as a ghost index at construction.
+       *
+       * This function does the same thing as operator().
+       */
+      Number &operator [] (const size_type global_index);
+
+      /**
+       * Read access to the data field specified by @p local_index. Locally
+       * owned indices can be accessed with indices
+       * <code>[0,local_size)</code>, and ghost indices with indices
+       * <code>[local_size,local_size+ n_ghost_entries]</code>.
+       *
+       * Performance: Direct array access (fast).
+       */
+      Number local_element (const size_type local_index) const;
+
+      /**
+       * Read and write access to the data field specified by @p local_index.
+       * Locally owned indices can be accessed with indices
+       * <code>[0,local_size)</code>, and ghost indices with indices
+       * <code>[local_size,local_size+n_ghosts]</code>.
+       *
+       * Performance: Direct array access (fast).
+       */
+      Number &local_element (const size_type local_index);
+
+      /**
+       * A collective get operation: instead of getting individual elements of
+       * a vector, this function allows to get a whole set of elements at
+       * once. The indices of the elements to be read are stated in the first
+       * argument, the corresponding values are returned in the second.
+       */
+      template <typename OtherNumber>
+      void extract_subvector_to (const std::vector<size_type> &indices,
+                                 std::vector<OtherNumber> &values) const;
+
+      /**
+       * Just as the above, but with pointers. Useful in minimizing copying of
+       * data around.
+       */
+      template <typename ForwardIterator, typename OutputIterator>
+      void extract_subvector_to (ForwardIterator          indices_begin,
+                                 const ForwardIterator    indices_end,
+                                 OutputIterator           values_begin) const;
+      /**
+       * Return whether the vector contains only elements with value zero.
+       * This is a collective operation. This function is expensive, because
+       * potentially all elements have to be checked.
+       */
+      bool all_zero () const;
+
+      /**
+       * Computes the mean value of all the entries in the vector.
+       */
+      Number mean_value () const;
+
+      /**
+       * $l_p$-norm of the vector. The pth root of the sum of the pth powers
+       * of the absolute values of the elements.
+       */
+      real_type lp_norm (const real_type p) const;
+      //@}
+
+      /**
+       * @name 6: Mixed stuff
+       */
+      //@{
+
+      /**
+       * Return a reference to the MPI communicator object in use with this
+       * vector.
+       */
+      const MPI_Comm &get_mpi_communicator () const;
+
+      /**
+       * Return the MPI partitioner that describes the parallel layout of the
+       * vector. This object can be used to initialize another vector with the
+       * respective reinit() call, for additional queries regarding the
+       * parallel communication, or the compatibility of partitioners.
+       */
+      std_cxx11::shared_ptr<const Utilities::MPI::Partitioner>
+      get_partitioner () const;
+
+      /**
+       * Checks whether the given partitioner is compatible with the
+       * partitioner used for this vector. Two partitioners are compatible if
+       * they have the same local size and the same ghost indices. They do not
+       * necessarily need to be the same data field of the shared pointer.
+       * This is a local operation only, i.e., if only some processors decide
+       * that the partitioning is not compatible, only these processors will
+       * return @p false, whereas the other processors will return @p true.
+       */
+      bool
+      partitioners_are_compatible (const Utilities::MPI::Partitioner &part) const;
+
+      /**
+       * Checks whether the given partitioner is compatible with the
+       * partitioner used for this vector. Two partitioners are compatible if
+       * they have the same local size and the same ghost indices. They do not
+       * necessarily need to be the same data field. As opposed to
+       * partitioners_are_compatible(), this method checks for compatibility
+       * among all processors and the method only returns @p true if the
+       * partitioner is the same on all processors.
+       *
+       * This method performs global communication, so make sure to use it
+       * only in a context where all processors call it the same number of
+       * times.
+       */
+      bool
+      partitioners_are_globally_compatible (const Utilities::MPI::Partitioner &part) const;
+      //@}
+
+      /**
+       * Attempt to perform an operation between two incompatible vector types.
+       *
+       * @ingroup Exceptions
+       */
+      DeclException0(ExcVectorTypeNotCompatible);
+
+      /**
+       * Exception
+       */
+      DeclException3 (ExcNonMatchingElements,
+                      double, double, unsigned int,
+                      << "Called compress(VectorOperation::insert), but"
+                      << " the element received from a remote processor, value "
+                      << std::setprecision(16) << arg1
+                      << ", does not match with the value "
+                      << std::setprecision(16) << arg2
+                      << " on the owner processor " << arg3);
+
+      /**
+       * Exception
+       */
+      DeclException4 (ExcAccessToNonLocalElement,
+                      size_type, size_type, size_type, size_type,
+                      << "You tried to access element " << arg1
+                      << " of a distributed vector, but this element is not "
+                      << "stored on the current processor. Note: The range of "
+                      << "locally owned elements is " << arg2 << " to "
+                      << arg3 << ", and there are " << arg4 << " ghost elements "
+                      << "that this vector can access.");
+
+    private:
+      /**
+       * Local part of all_zero().
+       */
+      bool all_zero_local () const;
+
+      /**
+       * Local part of the inner product of two vectors.
+       */
+      template <typename Number2>
+      Number inner_product_local (const Vector<Number2> &V) const;
+
+      /**
+       * Local part of norm_sqr().
+       */
+      real_type norm_sqr_local () const;
+
+      /**
+       * Local part of mean_value().
+       */
+      Number mean_value_local () const;
+
+      /**
+       * Local part of l1_norm().
+       */
+      real_type l1_norm_local () const;
+
+      /**
+       * Local part of lp_norm().
+       */
+      real_type lp_norm_local (const real_type p) const;
+
+      /**
+       * Local part of linfty_norm().
+       */
+      real_type linfty_norm_local () const;
+
+      /**
+       * Local part of the addition followed by an inner product of two
+       * vectors.
+       */
+      Number add_and_dot_local (const Number          a,
+                                const Vector<Number> &V,
+                                const Vector<Number> &W);
+
+      /**
+       * Shared pointer to store the parallel partitioning information. This
+       * information can be shared between several vectors that have the same
+       * partitioning.
+       */
+      std_cxx11::shared_ptr<const Utilities::MPI::Partitioner> partitioner;
+
+      /**
+       * The size that is currently allocated in the val array.
+       */
+      size_type allocated_size;
+
+      /**
+       * Pointer to the array of local elements of this vector.
+       */
+      Number         *val;
+
+      /**
+       * For parallel loops with TBB, this member variable stores the affinity
+       * information of loops.
+       */
+      mutable std_cxx11::shared_ptr<::dealii::parallel::internal::TBBPartitioner> thread_loop_partitioner;
+
+      /**
+       * Temporary storage that holds the data that is sent to this processor
+       * in @p compress() or sent from this processor in @p
+       * update_ghost_values.
+       */
+      mutable Number *import_data;
+
+      /**
+       * Stores whether the vector currently allows for reading ghost elements
+       * or not. Note that this is to ensure consistent ghost data and does
+       * not indicate whether the vector actually can store ghost elements. In
+       * particular, when assembling a vector we do not allow reading
+       * elements, only writing them.
+       */
+      mutable bool vector_is_ghosted;
+
+#ifdef DEAL_II_WITH_MPI
+      /**
+       * A vector that collects all requests from @p compress() operations.
+       * This class uses persistent MPI communicators, i.e., the communication
+       * channels are stored during successive calls to a given function. This
+       * reduces the overhead involved with setting up the MPI machinery, but
+       * it does not remove the need for a receive operation to be posted
+       * before the data can actually be sent.
+       */
+      std::vector<MPI_Request>   compress_requests;
+
+      /**
+       * A vector that collects all requests from @p update_ghost_values()
+       * operations. This class uses persistent MPI communicators.
+       */
+      mutable std::vector<MPI_Request>   update_ghost_values_requests;
+#endif
+
+      /**
+       * A lock that makes sure that the @p compress and @p
+       * update_ghost_values functions give reasonable results also when used
+       * with several threads.
+       */
+      mutable Threads::Mutex mutex;
+
+      /**
+       * A helper function that clears the compress_requests and
+       * update_ghost_values_requests field. Used in reinit functions.
+       */
+      void clear_mpi_requests ();
+
+      /**
+       * A helper function that is used to resize the val array.
+       */
+      void resize_val (const size_type new_allocated_size);
+
+      /*
+       * Make all other vector types friends.
+       */
+      template <typename Number2> friend class Vector;
+
+      /**
+       * Make BlockVector type friends.
+       */
+      template <typename Number2> friend class dealii::parallel::distributed::BlockVector;
+    };
+    /*@}*/
+
+
+    /*----------------------- Inline functions ----------------------------------*/
+
+#ifndef DOXYGEN
+
+
+    template <typename Number>
+    inline
+    bool
+    Vector<Number>::has_ghost_elements () const
+    {
+      return vector_is_ghosted;
+    }
+
+
+
+    template <typename Number>
+    inline
+    typename Vector<Number>::size_type
+    Vector<Number>::size () const
+    {
+      return partitioner->size();
+    }
+
+
+
+    template <typename Number>
+    inline
+    typename Vector<Number>::size_type
+    Vector<Number>::local_size () const
+    {
+      return partitioner->local_size();
+    }
+
+
+
+    template <typename Number>
+    inline
+    std::pair<typename Vector<Number>::size_type,
+        typename Vector<Number>::size_type>
+        Vector<Number>::local_range () const
+    {
+      return partitioner->local_range();
+    }
+
+
+
+    template <typename Number>
+    inline
+    bool
+    Vector<Number>::in_local_range
+    (const size_type global_index) const
+    {
+      return partitioner->in_local_range (global_index);
+    }
+
+
+
+    template <typename Number>
+    inline
+    IndexSet
+    Vector<Number>::locally_owned_elements() const
+    {
+      IndexSet is (size());
+
+      is.add_range (local_range().first, local_range().second);
+
+      return is;
+    }
+
+
+
+    template <typename Number>
+    inline
+    typename Vector<Number>::size_type
+    Vector<Number>::n_ghost_entries () const
+    {
+      return partitioner->n_ghost_indices();
+    }
+
+
+
+    template <typename Number>
+    inline
+    const IndexSet &
+    Vector<Number>::ghost_elements() const
+    {
+      return partitioner->ghost_indices();
+    }
+
+
+
+    template <typename Number>
+    inline
+    bool
+    Vector<Number>::is_ghost_entry (const size_type global_index) const
+    {
+      return partitioner->is_ghost_entry (global_index);
+    }
+
+
+
+    template <typename Number>
+    inline
+    typename Vector<Number>::iterator
+    Vector<Number>::begin ()
+    {
+      return val;
+    }
+
+
+
+    template <typename Number>
+    inline
+    typename Vector<Number>::const_iterator
+    Vector<Number>::begin () const
+    {
+      return val;
+    }
+
+
+
+    template <typename Number>
+    inline
+    typename Vector<Number>::iterator
+    Vector<Number>::end ()
+    {
+      return val+partitioner->local_size();
+    }
+
+
+
+    template <typename Number>
+    inline
+    typename Vector<Number>::const_iterator
+    Vector<Number>::end () const
+    {
+      return val+partitioner->local_size();
+    }
+
+
+
+    template <typename Number>
+    inline
+    Number
+    Vector<Number>::operator() (const size_type global_index) const
+    {
+      Assert (in_local_range (global_index) ||
+              partitioner->ghost_indices().is_element(global_index),
+              ExcAccessToNonLocalElement(global_index, local_range().first,
+                                         local_range().second,
+                                         partitioner->ghost_indices().n_elements()));
+      // do not allow reading a vector which is not in ghost mode
+      Assert (in_local_range (global_index) || vector_is_ghosted == true,
+              ExcMessage("You tried to read a ghost element of this vector, "
+                         "but it has not imported its ghost values."));
+      return val[partitioner->global_to_local(global_index)];
+    }
+
+
+
+    template <typename Number>
+    inline
+    Number &
+    Vector<Number>::operator() (const size_type global_index)
+    {
+      Assert (in_local_range (global_index) ||
+              partitioner->ghost_indices().is_element(global_index),
+              ExcAccessToNonLocalElement(global_index, local_range().first,
+                                         local_range().second,
+                                         partitioner->ghost_indices().n_elements()));
+      // we would like to prevent reading ghosts from a vector that does not
+      // have them imported, but this is not possible because we might be in a
+      // part of the code where the vector has enabled ghosts but is non-const
+      // (then, the compiler picks this method according to the C++ rule book
+      // even if a human would pick the const method when this subsequent use
+      // is just a read)
+      return val[partitioner->global_to_local (global_index)];
+    }
+
+
+
+    template <typename Number>
+    inline
+    Number
+    Vector<Number>::operator[] (const size_type global_index) const
+    {
+      return operator()(global_index);
+    }
+
+
+
+    template <typename Number>
+    inline
+    Number &
+    Vector<Number>::operator[] (const size_type global_index)
+    {
+      return operator()(global_index);
+    }
+
+
+
+    template <typename Number>
+    inline
+    Number
+    Vector<Number>::local_element (const size_type local_index) const
+    {
+      AssertIndexRange (local_index,
+                        partitioner->local_size()+
+                        partitioner->n_ghost_indices());
+      // do not allow reading a vector which is not in ghost mode
+      Assert (local_index < local_size() || vector_is_ghosted == true,
+              ExcMessage("You tried to read a ghost element of this vector, "
+                         "but it has not imported its ghost values."));
+      return val[local_index];
+    }
+
+
+
+    template <typename Number>
+    inline
+    Number &
+    Vector<Number>::local_element (const size_type local_index)
+    {
+      AssertIndexRange (local_index,
+                        partitioner->local_size()+
+                        partitioner->n_ghost_indices());
+      return val[local_index];
+    }
+
+
+
+    template <typename Number>
+    template <typename OtherNumber>
+    inline
+    void Vector<Number>::extract_subvector_to (const std::vector<size_type> &indices,
+                                               std::vector<OtherNumber> &values) const
+    {
+      for (size_type i = 0; i < indices.size(); ++i)
+        values[i] = operator()(indices[i]);
+    }
+
+
+
+    template <typename Number>
+    template <typename ForwardIterator, typename OutputIterator>
+    inline
+    void Vector<Number>::extract_subvector_to (ForwardIterator          indices_begin,
+                                               const ForwardIterator    indices_end,
+                                               OutputIterator           values_begin) const
+    {
+      while (indices_begin != indices_end)
+        {
+          *values_begin = operator()(*indices_begin);
+          indices_begin++;
+          values_begin++;
+        }
+    }
+
+
+
+    template <typename Number>
+    template <typename OtherNumber>
+    inline
+    void
+    Vector<Number>::add (const std::vector<size_type>    &indices,
+                         const ::dealii::Vector<OtherNumber> &values)
+    {
+      AssertDimension (indices.size(), values.size());
+      for (size_type i=0; i<indices.size(); ++i)
+        {
+          Assert (numbers::is_finite(values[i]),
+                  ExcMessage("The given value is not finite but either infinite or Not A Number (NaN)"));
+          this->operator()(indices[i]) += values(i);
+        }
+    }
+
+
+
+    template <typename Number>
+    inline
+    const MPI_Comm &
+    Vector<Number>::get_mpi_communicator() const
+    {
+      return partitioner->get_mpi_communicator();
+    }
+
+
+
+    template <typename Number>
+    inline
+    std_cxx11::shared_ptr<const Utilities::MPI::Partitioner>
+    Vector<Number>::get_partitioner () const
+    {
+      return partitioner;
+    }
+
+#endif
+
+  }
+}
+
+
+
+/**
+ * Global function @p swap which overloads the default implementation of the
+ * C++ standard library which uses a temporary object. The function simply
+ * exchanges the data of the two vectors.
+ *
+ * @relates Vector
+ * @author Katharina Kormann, Martin Kronbichler, 2011
+ */
+template <typename Number>
+inline
+void swap (LinearAlgebra::parallel::Vector<Number> &u,
+           LinearAlgebra::parallel::Vector<Number> &v)
+{
+  u.swap (v);
+}
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
similarity index 50%
rename from include/deal.II/lac/parallel_vector.templates.h
rename to include/deal.II/lac/la_parallel_vector.templates.h
index eca52c03d30e7ad63bc453d2f47d936d9bc16e59..2e6b7cf7ea531bcbb401d7a6cdd77f7490be8e42 100644 (file)
 //
 // ---------------------------------------------------------------------
 
-#ifndef dealii__parallel_vector_templates_h
-#define dealii__parallel_vector_templates_h
+#ifndef dealii__la_parallel_vector_templates_h
+#define dealii__la_parallel_vector_templates_h
 
 
 #include <deal.II/base/config.h>
-#include <deal.II/lac/parallel_vector.h>
-#include <deal.II/lac/vector_view.h>
+#include <deal.II/lac/la_parallel_vector.h>
+#include <deal.II/lac/vector_operations_internal.h>
+#include <deal.II/lac/read_write_vector.h>
 
-#include <deal.II/lac/petsc_parallel_vector.h>
-#include <deal.II/lac/trilinos_vector.h>
 
 DEAL_II_NAMESPACE_OPEN
 
 
-namespace parallel
+namespace LinearAlgebra
 {
-  namespace distributed
+  namespace parallel
   {
 
     template <typename Number>
@@ -70,6 +69,7 @@ namespace parallel
           val = 0;
           allocated_size = 0;
         }
+      thread_loop_partitioner.reset(new ::dealii::parallel::internal::TBBPartitioner());
     }
 
 
@@ -80,12 +80,10 @@ namespace parallel
                             const bool      omit_zeroing_entries)
     {
       clear_mpi_requests();
+
       // check whether we need to reallocate
       resize_val (size);
 
-      // reset vector view
-      vector_view.reinit (size, val);
-
       // delete previous content in import data
       if (import_data != 0)
         delete[] import_data;
@@ -122,11 +120,7 @@ namespace parallel
           const size_type new_allocated_size = partitioner->local_size() +
                                                partitioner->n_ghost_indices();
           resize_val (new_allocated_size);
-          vector_view.reinit (partitioner->local_size(), val);
         }
-      else
-        Assert (vector_view.size() == partitioner->local_size(),
-                ExcInternalError());
 
       if (omit_zeroing_entries == false)
         this->operator= (Number());
@@ -143,6 +137,8 @@ namespace parallel
         }
 
       vector_is_ghosted = false;
+
+      thread_loop_partitioner = v.thread_loop_partitioner;
     }
 
 
@@ -187,7 +183,6 @@ namespace parallel
       const size_type new_allocated_size = partitioner->local_size() +
                                            partitioner->n_ghost_indices();
       resize_val (new_allocated_size);
-      vector_view.reinit (partitioner->local_size(), val);
 
       // initialize to zero
       this->operator= (Number());
@@ -208,122 +203,215 @@ namespace parallel
 
 
 
-#ifdef DEAL_II_WITH_PETSC
+    template <typename Number>
+    Vector<Number>::Vector ()
+      :
+      partitioner (new Utilities::MPI::Partitioner()),
+      allocated_size (0),
+      val (0),
+      import_data (0)
+    {
+      reinit(0);
+    }
 
-    namespace internal
+
+
+    template <typename Number>
+    Vector<Number>::Vector (const Vector<Number> &v)
+      :
+      Subscriptor(),
+      allocated_size (0),
+      val (0),
+      import_data (0),
+      vector_is_ghosted (false)
     {
-      template <typename PETSC_Number, typename Number>
-      void copy_petsc_vector (const PETSC_Number *petsc_start_ptr,
-                              const PETSC_Number *petsc_end_ptr,
-                              Number *ptr)
-      {
-        std::copy(petsc_start_ptr, petsc_end_ptr, ptr);
-      }
+      reinit (v, true);
 
-      template <typename PETSC_Number, typename Number>
-      void copy_petsc_vector (const std::complex<PETSC_Number> *petsc_start_ptr,
-                              const std::complex<PETSC_Number> *petsc_end_ptr,
-                              std::complex<Number> *ptr)
-      {
-        std::copy(petsc_start_ptr, petsc_end_ptr, ptr);
-      }
+      thread_loop_partitioner = v.thread_loop_partitioner;
+      dealii::internal::Vector_copy<Number,Number> copier(v.val, val);
+      internal::parallel_for(copier, partitioner->local_size(),
+                             thread_loop_partitioner);
 
-      template <typename PETSC_Number, typename Number>
-      void copy_petsc_vector (const std::complex<PETSC_Number> * /*petsc_start_ptr*/,
-                              const std::complex<PETSC_Number> * /*petsc_end_ptr*/,
-                              Number * /*ptr*/)
-      {
-        AssertThrow(false, ExcMessage("Tried to copy complex -> real"));
-      }
+      zero_out_ghosts();
     }
 
+
+
     template <typename Number>
-    Vector<Number> &
-    Vector<Number>::operator = (const PETScWrappers::MPI::Vector &petsc_vec)
+    Vector<Number>::Vector (const IndexSet &local_range,
+                            const IndexSet &ghost_indices,
+                            const MPI_Comm  communicator)
+      :
+      allocated_size (0),
+      val (0),
+      import_data (0),
+      vector_is_ghosted (false)
     {
-      Assert(petsc_vec.locally_owned_elements() == locally_owned_elements(),
-             StandardExceptions::ExcInvalidState());
+      reinit (local_range, ghost_indices, communicator);
+    }
 
-      // get a representation of the vector and copy it
-      PetscScalar *start_ptr;
-      int ierr = VecGetArray (static_cast<const Vec &>(petsc_vec), &start_ptr);
-      AssertThrow (ierr == 0, ExcPETScError(ierr));
 
-      const size_type vec_size = local_size();
-      internal::copy_petsc_vector (start_ptr, start_ptr + vec_size, begin());
 
-      // restore the representation of the vector
-      ierr = VecRestoreArray (static_cast<const Vec &>(petsc_vec), &start_ptr);
-      AssertThrow (ierr == 0, ExcPETScError(ierr));
+    template <typename Number>
+    Vector<Number>::Vector (const IndexSet &local_range,
+                            const MPI_Comm  communicator)
+      :
+      allocated_size (0),
+      val (0),
+      import_data (0),
+      vector_is_ghosted (false)
+    {
+      reinit (local_range, communicator);
+    }
+
 
-      // spread ghost values between processes?
-      if (vector_is_ghosted || petsc_vec.has_ghost_elements())
-        update_ghost_values();
 
-      // return a pointer to this object per normal c++ operator overloading
-      // semantics
-      return *this;
+    template <typename Number>
+    Vector<Number>::Vector (const size_type size)
+      :
+      allocated_size (0),
+      val (0),
+      import_data (0),
+      vector_is_ghosted (false)
+    {
+      reinit (size, false);
     }
 
-#endif
 
 
+    template <typename Number>
+    Vector<Number>::
+    Vector (const std_cxx11::shared_ptr<const Utilities::MPI::Partitioner> &partitioner)
+      :
+      allocated_size (0),
+      val (0),
+      import_data (0),
+      vector_is_ghosted (false)
+    {
+      reinit (partitioner);
+    }
+
+
+
+    template <typename Number>
+    inline
+    Vector<Number>::~Vector ()
+    {
+      resize_val(0);
+
+      if (import_data != 0)
+        delete[] import_data;
+      import_data = 0;
+
+      clear_mpi_requests();
+    }
+
 
-#ifdef DEAL_II_WITH_TRILINOS
 
     template <typename Number>
+    inline
     Vector<Number> &
-    Vector<Number>::operator = (const TrilinosWrappers::MPI::Vector &trilinos_vec)
+    Vector<Number>::operator = (const Vector<Number> &c)
     {
-      if (trilinos_vec.has_ghost_elements() == false)
+#ifdef _MSC_VER
+      return this->operator=<Number>(c);
+#else
+      return this->template operator=<Number>(c);
+#endif
+    }
+
+
+
+    template <typename Number>
+    template <typename Number2>
+    inline
+    Vector<Number> &
+    Vector<Number>::operator = (const Vector<Number2> &c)
+    {
+      Assert (c.partitioner.get() != 0, ExcNotInitialized());
+
+      // we update ghost values whenever one of the input or output vector
+      // already held ghost values or when we import data from a vector with
+      // the same local range but different ghost layout
+      bool must_update_ghost_values = c.vector_is_ghosted;
+
+      // check whether the two vectors use the same parallel partitioner. if
+      // not, check if all local ranges are the same (that way, we can
+      // exchange data between different parallel layouts). One variant which
+      // is included here and necessary for compatibility with the other
+      // distributed vector classes (Trilinos, PETSc) is the case when vector
+      // c does not have any ghosts (constructed without ghost elements given)
+      // but the current vector does: In that case, we need to exchange data
+      // also when none of the two vector had updated its ghost values before.
+      if (partitioner.get() == 0)
+        reinit (c, true);
+      else if (partitioner.get() != c.partitioner.get())
         {
-          Assert(trilinos_vec.locally_owned_elements() == locally_owned_elements(),
-                 StandardExceptions::ExcInvalidState());
+          // local ranges are also the same if both partitioners are empty
+          // (even if they happen to define the empty range as [0,0) or [c,c)
+          // for some c!=0 in a different way).
+          int local_ranges_are_identical =
+            (local_range() == c.local_range() ||
+             (local_range().second == local_range().first &&
+              c.local_range().second == c.local_range().first));
+          if ((c.partitioner->n_mpi_processes() > 1 &&
+               Utilities::MPI::min(local_ranges_are_identical,
+                                   c.partitioner->get_communicator()) == 0)
+              ||
+              !local_ranges_are_identical)
+            reinit (c, true);
+          else
+            must_update_ghost_values |= vector_is_ghosted;
+
+          must_update_ghost_values |=
+            (c.partitioner->ghost_indices_initialized() == false &&
+             partitioner->ghost_indices_initialized() == true);
         }
       else
-        // ghosted trilinos vector must contain the local range of this vector
-        // which is contiguous
-        {
-          Assert((trilinos_vec.locally_owned_elements() & locally_owned_elements())
-                 == locally_owned_elements(),
-                 StandardExceptions::ExcInvalidState());
-        }
+        must_update_ghost_values |= vector_is_ghosted;
 
-      // create on trilinos data
-      const std::size_t start_index =
-        trilinos_vec.vector_partitioner().NumMyElements() > 0 ?
-        trilinos_vec.vector_partitioner().
-        LID(static_cast<TrilinosWrappers::types::int_type>(this->local_range().first)) : 0;
-      const VectorView<double> in_view (local_size(), trilinos_vec.begin()+start_index);
-      static_cast<dealii::Vector<Number>&>(vector_view) =
-        static_cast<const dealii::Vector<double>&>(in_view);
-
-      // spread ghost values between processes?
-      if (vector_is_ghosted || trilinos_vec.has_ghost_elements())
-        update_ghost_values();
+      thread_loop_partitioner = c.thread_loop_partitioner;
+      dealii::internal::Vector_copy<Number,Number2> copier(c.val, val);
+      internal::parallel_for(copier, partitioner->local_size(),
+                             thread_loop_partitioner);
 
-      // return a pointer to this object per normal c++ operator overloading
-      // semantics
+      if (must_update_ghost_values)
+        update_ghost_values();
+      else
+        zero_out_ghosts();
       return *this;
     }
 
-#endif
+
+
+    template <typename Number>
+    void
+    Vector<Number>::compress (::dealii::VectorOperation::values operation)
+    {
+      compress_start (0, operation);
+      compress_finish(operation);
+    }
 
 
 
     template <typename Number>
     void
-    Vector<Number>::copy_from (const Vector<Number> &c,
-                               const bool            call_update_ghost_values)
-    {
-      AssertDimension (local_range().first, c.local_range().first);
-      AssertDimension (local_range().second, c.local_range().second);
-      AssertDimension (vector_view.size(), c.vector_view.size());
-      vector_view = c.vector_view;
-      if (call_update_ghost_values == true)
-        update_ghost_values();
-      else
-        vector_is_ghosted = false;
+    Vector<Number>::update_ghost_values () const
+    {
+      update_ghost_values_start ();
+      update_ghost_values_finish ();
+    }
+
+
+
+    template <typename Number>
+    void
+    Vector<Number>::zero_out_ghosts ()
+    {
+      std::fill_n (&val[partitioner->local_size()],
+                   partitioner->n_ghost_indices(),
+                   Number());
+      vector_is_ghosted = false;
     }
 
 
@@ -396,7 +484,6 @@ namespace parallel
             }
           AssertDimension(current_index_start, part.n_import_indices());
 
-          Assert (part.local_size() == vector_view.size(), ExcInternalError());
           current_index_start = part.local_size();
           for (unsigned int i=0; i<n_ghost_targets; i++)
             {
@@ -539,8 +626,6 @@ namespace parallel
       // then actually send the data
       if (update_ghost_values_requests.size() == 0)
         {
-          Assert (part.local_size() == vector_view.size(),
-                  ExcInternalError());
           size_type current_index_start = part.local_size();
           update_ghost_values_requests.resize (n_import_targets+n_ghost_targets);
           for (unsigned int i=0; i<n_ghost_targets; i++)
@@ -633,6 +718,47 @@ namespace parallel
 
 
 
+    template <typename Number>
+    void
+    Vector<Number>::import(const ReadWriteVector<Number>                  &V,
+                           VectorOperation::values                         operation,
+                           std_cxx11::shared_ptr<const CommunicationPatternBase> communication_pattern)
+    {
+      // If no communication pattern is given, create one. Otherwise, use the
+      // given one.
+      std_cxx11::shared_ptr<const Utilities::MPI::Partitioner> comm_pattern;
+      if (communication_pattern.get() == NULL)
+        {
+          comm_pattern.reset(new Utilities::MPI::Partitioner(locally_owned_elements(),
+                                                             V.get_stored_elements(),
+                                                             get_mpi_communicator()));
+        }
+      else
+        {
+          comm_pattern =
+            std_cxx11::dynamic_pointer_cast<const Utilities::MPI::Partitioner> (communication_pattern);
+          AssertThrow(comm_pattern != NULL,
+                      ExcMessage(std::string("The communication pattern is not of type ") +
+                                 "Utilities::MPI::Partitioner."));
+        }
+      LinearAlgebra::parallel::Vector<Number> tmp_vector(comm_pattern);
+
+      // fill entries from ReadWriteVector into the distributed vector,
+      // including ghost entries. this is not really efficient right now
+      // because indices are translated twice, once for
+      const IndexSet &v_stored = V.get_stored_elements();
+      for (size_type i=0; i<v_stored.n_elements(); ++i)
+        tmp_vector(v_stored.nth_index_in_set(i)) = V.local_element(i);
+
+      tmp_vector.compress(operation);
+
+      dealii::internal::Vector_copy<Number,Number> copier(tmp_vector.val, val);
+      internal::parallel_for(copier, partitioner->local_size(),
+                             thread_loop_partitioner);
+    }
+
+
+
     template <typename Number>
     void
     Vector<Number>::swap (Vector<Number> &v)
@@ -672,15 +798,575 @@ namespace parallel
 #endif
 
       std::swap (partitioner,       v.partitioner);
+      std::swap (thread_loop_partitioner, v.thread_loop_partitioner);
       std::swap (allocated_size,    v.allocated_size);
       std::swap (val,               v.val);
       std::swap (import_data,       v.import_data);
       std::swap (vector_is_ghosted, v.vector_is_ghosted);
+    }
+
+
+
+    template <typename Number>
+    Vector<Number> &
+    Vector<Number>::operator = (const Number s)
+    {
+      internal::Vector_set<Number> setter(s, val);
+
+      internal::parallel_for(setter, partitioner->local_size(),
+                             thread_loop_partitioner);
+
+      // if we call Vector::operator=0, we want to zero out all the entries
+      // plus ghosts.
+      if (s==Number())
+        zero_out_ghosts();
+
+      return *this;
+    }
+
+
+
+    template <typename Number>
+    Vector<Number> &
+    Vector<Number>::operator += (const VectorSpaceVector<Number> &vv)
+    {
+      // Downcast. Throws an exception if invalid.
+      Assert(dynamic_cast<const Vector<Number> *>(&vv)!=NULL,
+             ExcVectorTypeNotCompatible());
+      const Vector<Number> &v = dynamic_cast<const Vector<Number> &>(vv);
+
+      AssertDimension (local_size(), v.local_size());
+
+      internal::Vectorization_add_v<Number> vector_add(val, v.val);
+      internal::parallel_for(vector_add, partitioner->local_size(),
+                             thread_loop_partitioner);
+
+      if (vector_is_ghosted)
+        update_ghost_values();
+
+      return *this;
+    }
+
+
+
+    template <typename Number>
+    Vector<Number> &
+    Vector<Number>::operator -= (const VectorSpaceVector<Number> &vv)
+    {
+      // Downcast. Throws an exception if invalid.
+      Assert(dynamic_cast<const Vector<Number> *>(&vv)!=NULL,
+             ExcVectorTypeNotCompatible());
+      const Vector<Number> &v = dynamic_cast<const Vector<Number> &>(vv);
+
+      AssertDimension (local_size(), v.local_size());
+
+      internal::Vectorization_subtract_v<Number> vector_subtract(val, v.val);
+      internal::parallel_for(vector_subtract, partitioner->local_size(),
+                             thread_loop_partitioner);
+
+      if (vector_is_ghosted)
+        update_ghost_values();
+
+      return *this;
+    }
+
+
+
+    template <typename Number>
+    void
+    Vector<Number>::add (const Number a)
+    {
+      AssertIsFinite(a);
+
+      internal::Vectorization_add_factor<Number> vector_add(val, a);
+      internal::parallel_for(vector_add, partitioner->local_size(),
+                             thread_loop_partitioner);
+
+      if (vector_is_ghosted)
+        update_ghost_values();
+    }
+
+
+
+    template <typename Number>
+    void
+    Vector<Number>::add (const Number a,
+                         const VectorSpaceVector<Number> &vv)
+    {
+      // Downcast. Throws an exception if invalid.
+      Assert(dynamic_cast<const Vector<Number> *>(&vv)!=NULL,
+             ExcVectorTypeNotCompatible());
+      const Vector<Number> &v = dynamic_cast<const Vector<Number> &>(vv);
+
+      AssertIsFinite(a);
+      AssertDimension (local_size(), v.local_size());
+
+      internal::Vectorization_add_av<Number> vector_add(val, v.val, a);
+      internal::parallel_for(vector_add, partitioner->local_size(),
+                             thread_loop_partitioner);
+
+      if (vector_is_ghosted)
+        update_ghost_values();
+    }
+
+
+
+    template <typename Number>
+    void
+    Vector<Number>::add (const Number a,
+                         const VectorSpaceVector<Number> &vv,
+                         const Number b,
+                         const VectorSpaceVector<Number> &ww)
+    {
+      // Downcast. Throws an exception if invalid.
+      Assert(dynamic_cast<const Vector<Number> *>(&vv)!=NULL,
+             ExcVectorTypeNotCompatible());
+      const Vector<Number> &v = dynamic_cast<const Vector<Number> &>(vv);
+      Assert(dynamic_cast<const Vector<Number> *>(&ww)!=NULL,
+             ExcVectorTypeNotCompatible());
+      const Vector<Number> &w = dynamic_cast<const Vector<Number> &>(ww);
+
+      AssertIsFinite(a);
+      AssertIsFinite(b);
+
+      AssertDimension (local_size(), v.local_size());
+      AssertDimension (local_size(), w.local_size());
+
+      internal::Vectorization_add_avpbw<Number> vector_add(val, v.val, w.val, a, b);
+      internal::parallel_for(vector_add, partitioner->local_size(),
+                             thread_loop_partitioner);
+
+      if (vector_is_ghosted)
+        update_ghost_values();
+    }
+
+
+
+    template <typename Number>
+    void
+    Vector<Number>::sadd (const Number x,
+                          const Vector<Number> &v)
+    {
+      AssertIsFinite(x);
+      AssertDimension (local_size(), v.local_size());
+
+      internal::Vectorization_sadd_xv<Number> vector_sadd(val, v.val, x);
+      internal::parallel_for(vector_sadd, partitioner->local_size(),
+                             thread_loop_partitioner);
+
+      if (vector_is_ghosted)
+        update_ghost_values();
+    }
+
+
+
+    template <typename Number>
+    void
+    Vector<Number>::sadd (const Number x,
+                          const Number a,
+                          const VectorSpaceVector<Number> &vv)
+    {
+      // Downcast. Throws an exception if invalid.
+      Assert(dynamic_cast<const Vector<Number> *>(&vv)!=NULL,
+             ExcVectorTypeNotCompatible());
+      const Vector<Number> &v = dynamic_cast<const Vector<Number> &>(vv);
+
+      AssertIsFinite(x);
+      AssertIsFinite(a);
+      AssertDimension (local_size(), v.local_size());
+
+      internal::Vectorization_sadd_xav<Number> vector_sadd(val, v.val, a, x);
+      internal::parallel_for(vector_sadd, partitioner->local_size(),
+                             thread_loop_partitioner);
+
+      if (vector_is_ghosted)
+        update_ghost_values();
+    }
 
-      // vector view cannot be swapped so reset it manually (without touching
-      // the vector elements)
-      vector_view.reinit (partitioner->local_size(), val);
-      v.vector_view.reinit (v.partitioner->local_size(), v.val);
+
+
+    template <typename Number>
+    void
+    Vector<Number>::sadd (const Number x,
+                          const Number a,
+                          const Vector<Number> &v,
+                          const Number b,
+                          const Vector<Number> &w)
+    {
+      AssertIsFinite(x);
+      AssertIsFinite(a);
+      AssertIsFinite(b);
+
+      AssertDimension (local_size(), v.local_size());
+      AssertDimension (local_size(), w.local_size());
+
+      internal::Vectorization_sadd_xavbw<Number> vector_sadd(val, v.val, w.val,
+                                                             x, a, b);
+      internal::parallel_for(vector_sadd, partitioner->local_size(),
+                             thread_loop_partitioner);
+
+      if (vector_is_ghosted)
+        update_ghost_values();
+    }
+
+
+
+    template <typename Number>
+    Vector<Number> &
+    Vector<Number>::operator *= (const Number factor)
+    {
+      AssertIsFinite(factor);
+      internal::Vectorization_multiply_factor<Number> vector_multiply(val, factor);
+
+      internal::parallel_for(vector_multiply, partitioner->local_size(),
+                             thread_loop_partitioner);
+
+      if (vector_is_ghosted)
+        update_ghost_values();
+
+      return *this;
+    }
+
+
+
+    template <typename Number>
+    Vector<Number> &
+    Vector<Number>::operator /= (const Number factor)
+    {
+      operator *= (1./factor);
+      return *this;
+    }
+
+
+
+    template <typename Number>
+    void
+    Vector<Number>::scale (const VectorSpaceVector<Number> &vv)
+    {
+      // Downcast. Throws an exception if invalid.
+      Assert(dynamic_cast<const Vector<Number> *>(&vv)!=NULL,
+             ExcVectorTypeNotCompatible());
+      const Vector<Number> &v = dynamic_cast<const Vector<Number> &>(vv);
+
+      AssertDimension (local_size(), v.local_size());
+
+      internal::Vectorization_scale<Number> vector_scale(val, v.val);
+      internal::parallel_for(vector_scale, partitioner->local_size(),
+                             thread_loop_partitioner);
+
+      if (vector_is_ghosted)
+        update_ghost_values();
+    }
+
+
+
+    template <typename Number>
+    void
+    Vector<Number>::equ (const Number a,
+                         const VectorSpaceVector<Number> &vv)
+    {
+      // Downcast. Throws an exception if invalid.
+      Assert(dynamic_cast<const Vector<Number> *>(&vv)!=NULL,
+             ExcVectorTypeNotCompatible());
+      const Vector<Number> &v = dynamic_cast<const Vector<Number> &>(vv);
+
+      AssertIsFinite(a);
+      AssertDimension (local_size(), v.local_size());
+
+      internal::Vectorization_equ_au<Number> vector_equ(val, v.val, a);
+      internal::parallel_for(vector_equ, partitioner->local_size(),
+                             thread_loop_partitioner);
+
+      if (vector_is_ghosted)
+        update_ghost_values();
+    }
+
+
+
+    template <typename Number>
+    void
+    Vector<Number>::equ (const Number a,
+                         const Vector<Number> &v,
+                         const Number b,
+                         const Vector<Number> &w)
+    {
+      AssertIsFinite(a);
+      AssertIsFinite(b);
+
+      AssertDimension (local_size(), v.local_size());
+      AssertDimension (local_size(), w.local_size());
+
+      internal::Vectorization_equ_aubv<Number> vector_equ(val, v.val, w.val, a, b);
+      internal::parallel_for(vector_equ, partitioner->local_size(),
+                             thread_loop_partitioner);
+
+      if (vector_is_ghosted)
+        update_ghost_values();
+    }
+
+
+
+    template <typename Number>
+    bool
+    Vector<Number>::all_zero_local () const
+    {
+      const size_type local_size = partitioner->local_size();
+      for (size_type i=0; i<local_size; ++i)
+        if (val[i] != Number(0))
+          return false;
+      return true;
+    }
+
+
+
+    template <typename Number>
+    bool
+    Vector<Number>::all_zero () const
+    {
+      // use int instead of bool. in order to make global reduction operations
+      // work also when MPI_Init was not called, only call MPI_Allreduce
+      // commands when there is more than one processor (note that reinit()
+      // functions handle this case correctly through the job_supports_mpi()
+      // query). this is the same in all the reduce functions below
+      int local_result = -static_cast<int>(all_zero_local());
+      if (partitioner->n_mpi_processes() > 1)
+        return -Utilities::MPI::max(local_result,
+                                    partitioner->get_communicator());
+      else
+        return -local_result;
+    }
+
+
+
+    template <typename Number>
+    template <typename Number2>
+    Number
+    Vector<Number>::inner_product_local(const Vector<Number2> &v) const
+    {
+      if (PointerComparison::equal (this, &v))
+        return norm_sqr_local();
+
+      AssertDimension (partitioner->local_size(), v.partitioner->local_size());
+
+      Number sum;
+      internal::Dot<Number,Number2> dot(val, v.val);
+      internal::parallel_reduce (dot, partitioner->local_size(), sum,
+                                 thread_loop_partitioner);
+      AssertIsFinite(sum);
+
+      return sum;
+    }
+
+
+
+    template <typename Number>
+    Number
+    Vector<Number>::operator * (const VectorSpaceVector<Number> &vv) const
+    {
+      // Downcast. Throws an exception if invalid.
+      Assert(dynamic_cast<const Vector<Number> *>(&vv)!=NULL,
+             ExcVectorTypeNotCompatible());
+      const Vector<Number> &v = dynamic_cast<const Vector<Number> &>(vv);
+
+      Number local_result = inner_product_local(v);
+      if (partitioner->n_mpi_processes() > 1)
+        return Utilities::MPI::sum (local_result,
+                                    partitioner->get_communicator());
+      else
+        return local_result;
+    }
+
+
+
+    template <typename Number>
+    typename Vector<Number>::real_type
+    Vector<Number>::norm_sqr_local () const
+    {
+      real_type sum;
+      internal::Norm2<Number,real_type> norm2(val);
+      internal::parallel_reduce (norm2, partitioner->local_size(), sum,
+                                 thread_loop_partitioner);
+      AssertIsFinite(sum);
+
+      return sum;
+    }
+
+
+
+    template <typename Number>
+    Number
+    Vector<Number>::mean_value_local () const
+    {
+      Assert(size() != 0, ExcEmptyObject());
+
+      if (partitioner->local_size() == 0)
+        return Number();
+
+      Number sum;
+      internal::MeanValue<Number> mean(val);
+      internal::parallel_reduce (mean, partitioner->local_size(), sum,
+                                 thread_loop_partitioner);
+
+      return sum / real_type(partitioner->local_size());
+    }
+
+
+
+    template <typename Number>
+    Number
+    Vector<Number>::mean_value () const
+    {
+      Number local_result = mean_value_local();
+      if (partitioner->n_mpi_processes() > 1)
+        return Utilities::MPI::sum (local_result *
+                                    (real_type)partitioner->local_size(),
+                                    partitioner->get_communicator())
+               /(real_type)partitioner->size();
+      else
+        return local_result;
+    }
+
+
+
+    template <typename Number>
+    typename Vector<Number>::real_type
+    Vector<Number>::l1_norm_local () const
+    {
+      real_type sum;
+      internal::Norm1<Number, real_type> norm1(val);
+      internal::parallel_reduce (norm1, partitioner->local_size(), sum,
+                                 thread_loop_partitioner);
+
+      return sum;
+    }
+
+
+
+    template <typename Number>
+    typename Vector<Number>::real_type
+    Vector<Number>::l1_norm () const
+    {
+      real_type local_result = l1_norm_local();
+      if (partitioner->n_mpi_processes() > 1)
+        return Utilities::MPI::sum(local_result,
+                                   partitioner->get_communicator());
+      else
+        return local_result;
+    }
+
+
+
+    template <typename Number>
+    typename Vector<Number>::real_type
+    Vector<Number>::l2_norm () const
+    {
+      real_type local_result = norm_sqr_local();
+      if (partitioner->n_mpi_processes() > 1)
+        return std::sqrt(Utilities::MPI::sum(local_result,
+                                             partitioner->get_communicator()));
+      else
+        return std::sqrt(local_result);
+    }
+
+
+
+    template <typename Number>
+    typename Vector<Number>::real_type
+    Vector<Number>::lp_norm_local (const real_type p) const
+    {
+      real_type sum;
+      internal::NormP<Number, real_type> normp(val, p);
+      internal::parallel_reduce (normp, partitioner->local_size(), sum,
+                                 thread_loop_partitioner);
+      return std::pow(sum, 1./p);
+    }
+
+
+
+    template <typename Number>
+    typename Vector<Number>::real_type
+    Vector<Number>::lp_norm (const real_type p) const
+    {
+      const real_type local_result = lp_norm_local(p);
+      if (partitioner->n_mpi_processes() > 1)
+        return std::pow (Utilities::MPI::sum(std::pow(local_result,p),
+                                             partitioner->get_communicator()),
+                         static_cast<real_type>(1.0/p));
+      else
+        return local_result;
+    }
+
+
+
+    template <typename Number>
+    typename Vector<Number>::real_type
+    Vector<Number>::linfty_norm_local () const
+    {
+      real_type max = 0.;
+
+      const size_type local_size = partitioner->local_size();
+      for (size_type i=0; i<local_size; ++i)
+        max = std::max (numbers::NumberTraits<Number>::abs(val[i]), max);
+
+      return max;
+    }
+
+
+
+    template <typename Number>
+    inline
+    typename Vector<Number>::real_type
+    Vector<Number>::linfty_norm () const
+    {
+      const real_type local_result = linfty_norm_local();
+      if (partitioner->n_mpi_processes() > 1)
+        return Utilities::MPI::max (local_result,
+                                    partitioner->get_communicator());
+      else
+        return local_result;
+    }
+
+
+
+    template <typename Number>
+    Number
+    Vector<Number>::add_and_dot_local(const Number          a,
+                                      const Vector<Number> &v,
+                                      const Vector<Number> &w)
+    {
+      const size_type vec_size = partitioner->local_size();
+      AssertDimension (vec_size, v.local_size());
+      AssertDimension (vec_size, w.local_size());
+
+      Number sum;
+      internal::AddAndDot<Number> adder(this->val, v.val, w.val, a);
+      internal::parallel_reduce (adder, vec_size, sum, thread_loop_partitioner);
+      AssertIsFinite(sum);
+      return sum;
+    }
+
+
+
+    template <typename Number>
+    Number
+    Vector<Number>::add_and_dot (const Number                     a,
+                                 const VectorSpaceVector<Number> &vv,
+                                 const VectorSpaceVector<Number> &ww)
+    {
+      // Downcast. Throws an exception if invalid.
+      Assert(dynamic_cast<const Vector<Number> *>(&vv)!=NULL,
+             ExcVectorTypeNotCompatible());
+      const Vector<Number> &v = dynamic_cast<const Vector<Number> &>(vv);
+      Assert(dynamic_cast<const Vector<Number> *>(&ww)!=NULL,
+             ExcVectorTypeNotCompatible());
+      const Vector<Number> &w = dynamic_cast<const Vector<Number> &>(ww);
+
+      Number local_result = add_and_dot_local(a, v, w);
+      if (partitioner->n_mpi_processes() > 1)
+        return Utilities::MPI::sum (local_result,
+                                    partitioner->get_communicator());
+      else
+        return local_result;
     }
 
 
index 7d472ae2df421aab31b00757649053fa9e179655..8cb7113f4258182b635c8e594ef0f55a6bffbf46 100644 (file)
 #define dealii__parallel_vector_h
 
 #include <deal.II/base/config.h>
-#include <deal.II/base/index_set.h>
-#include <deal.II/base/template_constraints.h>
-#include <deal.II/base/types.h>
-#include <deal.II/base/utilities.h>
-#include <deal.II/base/memory_consumption.h>
-#include <deal.II/base/partitioner.h>
-#include <deal.II/base/thread_management.h>
-#include <deal.II/lac/vector_view.h>
+#include <deal.II/lac/la_parallel_vector.h>
 
 #include <cstring>
 #include <iomanip>
 
 DEAL_II_NAMESPACE_OPEN
 
-#ifdef DEAL_II_WITH_PETSC
-namespace PETScWrappers
-{
-  namespace MPI
-  {
-    class Vector;
-  }
-}
-#endif
-
-#ifdef DEAL_II_WITH_TRILINOS
-namespace TrilinosWrappers
-{
-  namespace MPI
-  {
-    class Vector;
-  }
-}
-#endif
-
 
 namespace parallel
 {
@@ -166,2272 +139,10 @@ namespace parallel
      *
      * @author Katharina Kormann, Martin Kronbichler, 2010, 2011
      */
-    template <typename Number>
-    class Vector : public Subscriptor
-    {
-    public:
-      /**
-       * Declare standard types used in all containers. These types parallel
-       * those in the <tt>C++</tt> standard libraries <tt>vector<...></tt>
-       * class.
-       */
-      typedef Number                                            value_type;
-      typedef value_type                                       *pointer;
-      typedef const value_type                                 *const_pointer;
-      typedef value_type                                       *iterator;
-      typedef const value_type                                 *const_iterator;
-      typedef value_type                                       &reference;
-      typedef const value_type                                 &const_reference;
-      typedef types::global_dof_index                           size_type;
-      typedef typename numbers::NumberTraits<Number>::real_type real_type;
-
-      /**
-       * A variable that indicates whether this vector supports distributed
-       * data storage. If true, then this vector also needs an appropriate
-       * compress() function that allows communicating recent set or add
-       * operations to individual elements to be communicated to other
-       * processors.
-       *
-       * For the current class, the variable equals true, since it does
-       * support parallel data storage.
-       */
-      static const bool supports_distributed_data = true;
-
-      /**
-       * @name 1: Basic Object-handling
-       */
-      //@{
-      /**
-       * Empty constructor.
-       */
-      Vector ();
-
-      /**
-       * Copy constructor. Uses the parallel partitioning of @p in_vector.
-       */
-      Vector (const Vector<Number> &in_vector);
-
-      /**
-       * Constructs a parallel vector of the given global size without any
-       * actual parallel distribution.
-       */
-      Vector (const size_type size);
-
-      /**
-       * Constructs a parallel vector. The local range is specified by @p
-       * locally_owned_set (note that this must be a contiguous interval,
-       * multiple intervals are not possible). The IndexSet @p ghost_indices
-       * specifies ghost indices, i.e., indices which one might need to read
-       * data from or accumulate data from. It is allowed that the set of
-       * ghost indices also contains the local range, but it does not need to.
-       *
-       * This function involves global communication, so it should only be
-       * called once for a given layout. Use the constructor with
-       * Vector<Number> argument to create additional vectors with the same
-       * parallel layout.
-       *
-       * @see
-       * @ref GlossGhostedVector "vectors with ghost elements"
-       */
-      Vector (const IndexSet &local_range,
-              const IndexSet &ghost_indices,
-              const MPI_Comm  communicator);
-
-      /**
-       * Same constructor as above but without any ghost indices.
-       */
-      Vector (const IndexSet &local_range,
-              const MPI_Comm  communicator);
-
-      /**
-       * Create the vector based on the parallel partitioning described in @p
-       * partitioner. The input argument is a shared pointer, which store the
-       * partitioner data only once and share it between several vectors with
-       * the same layout.
-       */
-      Vector (const std_cxx11::shared_ptr<const Utilities::MPI::Partitioner> &partitioner);
-
-      /**
-       * Destructor.
-       */
-      ~Vector ();
-
-      /**
-       * Sets the global size of the vector to @p size without any actual
-       * parallel distribution.
-       */
-      void reinit (const size_type size,
-                   const bool      omit_zeroing_entries = false);
-
-      /**
-       * Uses the parallel layout of the input vector @p in_vector and
-       * allocates memory for this vector. Recommended initialization function
-       * when several vectors with the same layout should be created.
-       *
-       * 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 Vector<Number2> &in_vector,
-                  const bool             omit_zeroing_entries = false);
-
-      /**
-       * Initialize the vector. The local range is specified by @p
-       * locally_owned_set (note that this must be a contiguous interval,
-       * multiple intervals are not possible). The IndexSet @p ghost_indices
-       * specifies ghost indices, i.e., indices which one might need to read
-       * data from or accumulate data from. It is allowed that the set of
-       * ghost indices also contains the local range, but it does not need to.
-       *
-       * This function involves global communication, so it should only be
-       * called once for a given layout. Use the @p reinit function with
-       * Vector<Number> argument to create additional vectors with the same
-       * parallel layout.
-       *
-       * @see
-       * @ref GlossGhostedVector "vectors with ghost elements"
-       */
-      void reinit (const IndexSet &local_range,
-                   const IndexSet &ghost_indices,
-                   const MPI_Comm  communicator);
-
-      /**
-       * Same as above, but without ghost entries.
-       */
-      void reinit (const IndexSet &local_range,
-                   const MPI_Comm  communicator);
-
-      /**
-       * Initialize the vector given to the parallel partitioning described in
-       * @p partitioner. The input argument is a shared pointer, which store
-       * the partitioner data only once and share it between several vectors
-       * with the same layout.
-       */
-      void reinit (const std_cxx11::shared_ptr<const Utilities::MPI::Partitioner> &partitioner);
-
-      /**
-       * Swap the contents of this vector and the other vector @p v. One could
-       * do this operation with a temporary variable and copying over the data
-       * elements, but this function is significantly more efficient since it
-       * only swaps the pointers to the data of the two vectors and therefore
-       * does not need to allocate temporary storage and move data around.
-       *
-       * This function is analog to the the @p swap function of all C++
-       * standard containers. Also, there is a global function
-       * <tt>swap(u,v)</tt> that simply calls <tt>u.swap(v)</tt>, again in
-       * analogy to standard functions.
-       *
-       * This function is virtual in order to allow for derived classes to
-       * handle memory separately.
-       */
-      void swap (Vector<Number> &v);
-
-      /**
-       * Assigns the vector to the parallel partitioning of the input vector
-       * @p in_vector, and copies all the data.
-       *
-       * If one of the input vector or the calling vector (to the left of the
-       * assignment operator) had ghost elements set before this operation,
-       * the calling vector will have ghost values set. Otherwise, it will be
-       * in write mode. If the input vector does not have any ghost elements
-       * at all, the vector will also update its ghost values in analogy to
-       * the respective setting the Trilinos and PETSc vectors.
-       */
-      Vector<Number> &
-      operator = (const Vector<Number> &in_vector);
-
-      /**
-       * Assigns the vector to the parallel partitioning of the input vector
-       * @p in_vector, and copies all the data.
-       *
-       * If one of the input vector or the calling vector (to the left of the
-       * assignment operator) had ghost elements set before this operation,
-       * the calling vector will have ghost values set. Otherwise, it will be
-       * in write mode. If the input vector does not have any ghost elements
-       * at all, the vector will also update its ghost values in analogy to
-       * the respective setting the Trilinos and PETSc vectors.
-       */
-      template <typename Number2>
-      Vector<Number> &
-      operator = (const Vector<Number2> &in_vector);
-
-#ifdef DEAL_II_WITH_PETSC
-      /**
-       * Copy the content of a PETSc vector into the calling vector. This
-       * function assumes that the vectors layouts have already been
-       * initialized to match.
-       *
-       * This operator is only available if deal.II was configured with PETSc.
-       */
-      Vector<Number> &
-      operator = (const PETScWrappers::MPI::Vector &petsc_vec);
-#endif
-
-#ifdef DEAL_II_WITH_TRILINOS
-      /**
-       * Copy the content of a Trilinos vector into the calling vector. This
-       * function assumes that the vectors layouts have already been
-       * initialized to match.
-       *
-       * This operator is only available if deal.II was configured with
-       * Trilinos.
-       */
-      Vector<Number> &
-      operator = (const TrilinosWrappers::MPI::Vector &trilinos_vec);
-#endif
-
-      /**
-       * This method copies the local range from another vector with the same
-       * local range, but possibly different layout of ghost indices.
-       *
-       * This function is deprecated.
-       */
-      void copy_from (const Vector<Number> &in_vector,
-                      const bool            call_update_ghost_values = false) DEAL_II_DEPRECATED;
-
-      /**
-       * Sets all elements of the vector to the scalar @p s. If the scalar is
-       * zero, also ghost elements are set to zero, otherwise they remain
-       * unchanged.
-       */
-      Vector<Number> &operator = (const Number s);
-
-      /**
-       * This function copies the data that has accumulated in the data buffer
-       * for ghost indices to the owning processor. For the meaning of the
-       * argument @p operation, see the entry on
-       * @ref GlossCompress "Compressing distributed vectors and matrices"
-       * in the glossary.
-       *
-       * There are two variants for this function. If called with argument @p
-       * VectorOperation::add adds all the data accumulated in ghost elements
-       * to the respective elements on the owning processor and clears the
-       * ghost array afterwards. If called with argument @p
-       * VectorOperation::insert, a set operation is performed. Since setting
-       * elements in a vector with ghost elements is ambiguous (as one can set
-       * both the element on the ghost site as well as the owning site), this
-       * operation makes the assumption that all data is set correctly on the
-       * owning processor. Upon call of compress(VectorOperation::insert), all
-       * ghost entries are thus simply zeroed out (using zero_ghost_values()).
-       * In debug mode, a check is performed for whether the data set is
-       * actually consistent between processors, i.e., whenever a non-zero
-       * ghost element is found, it is compared to the value on the owning
-       * processor and an exception is thrown if these elements do not agree.
-       */
-      void compress (::dealii::VectorOperation::values operation);
-
-      /**
-       * Fills the data field for ghost indices with the values stored in the
-       * respective positions of the owning processor. This function is needed
-       * before reading from ghosts. The function is @p const even though
-       * ghost data is changed. This is needed to allow functions with a @p
-       * const vector to perform the data exchange without creating
-       * temporaries.
-       *
-       * After calling this method, write access to ghost elements of the
-       * vector is forbidden and an exception is thrown. Only read access to
-       * ghost elements is allowed in this state. Note that all subsequent
-       * operations on this vector, like global vector addition, etc., will
-       * also update the ghost values by a call to this method after the
-       * operation. However, global reduction operations like norms or the
-       * inner product will always ignore ghost elements in order to avoid
-       * counting the ghost data more than once. To allow writing to ghost
-       * elements again, call zero_out_ghosts().
-       *
-       * @see
-       * @ref GlossGhostedVector "vectors with ghost elements"
-       */
-      void update_ghost_values () const;
-
-      /**
-       * Initiates communication for the @p compress() function with non-
-       * blocking communication. This function does not wait for the transfer
-       * to finish, in order to allow for other computations during the time
-       * it takes until all data arrives.
-       *
-       * Before the data is actually exchanged, the function must be followed
-       * by a call to @p compress_finish().
-       *
-       * In case this function is called for more than one vector before @p
-       * compress_finish() is invoked, it is mandatory to specify a unique
-       * communication channel to each such call, in order to avoid several
-       * messages with the same ID that will corrupt this operation.
-       */
-      void compress_start (const unsigned int communication_channel = 0,
-                           ::dealii::VectorOperation::values operation = VectorOperation::add);
-
-      /**
-       * For all requests that have been initiated in compress_start, wait for
-       * the communication to finish. Once it is finished, add or set the data
-       * (depending on the flag operation) to the respective positions in the
-       * owning processor, and clear the contents in the ghost data fields.
-       * The meaning of this argument is the same as in compress().
-       *
-       * This function should be called exactly once per vector after calling
-       * compress_start, otherwise the result is undefined. In particular, it
-       * is not well-defined to call compress_start on the same vector again
-       * before compress_finished has been called. However, there is no
-       * warning to prevent this situation.
-       *
-       * Must follow a call to the @p compress_start function.
-       */
-      void compress_finish (::dealii::VectorOperation::values operation);
-
-      /**
-       * Initiates communication for the @p update_ghost_values() function
-       * with non-blocking communication. This function does not wait for the
-       * transfer to finish, in order to allow for other computations during
-       * the time it takes until all data arrives.
-       *
-       * Before the data is actually exchanged, the function must be followed
-       * by a call to @p update_ghost_values_finish().
-       *
-       * In case this function is called for more than one vector before @p
-       * update_ghost_values_finish() is invoked, it is mandatory to specify a
-       * unique communication channel to each such call, in order to avoid
-       * several messages with the same ID that will corrupt this operation.
-       */
-      void update_ghost_values_start (const unsigned int communication_channel = 0) const;
-
-
-      /**
-       * For all requests that have been started in update_ghost_values_start,
-       * wait for the communication to finish.
-       *
-       * Must follow a call to the @p update_ghost_values_start function
-       * before reading data from ghost indices.
-       */
-      void update_ghost_values_finish () const;
-
-      /**
-       * This method zeros the entries on ghost dofs, but does not touch
-       * locally owned DoFs.
-       *
-       * After calling this method, read access to ghost elements of the
-       * vector is forbidden and an exception is thrown. Only write access to
-       * ghost elements is allowed in this state.
-       */
-      void zero_out_ghosts ();
-
-      /**
-       * Returns whether the vector currently is in a state where ghost values
-       * can be read or not. This is the same functionality as other parallel
-       * vectors have. If this method returns false, this only means that
-       * read-access to ghost elements is prohibited whereas write access is
-       * still possible (to those entries specified as ghosts during
-       * initialization), not that there are no ghost elements at all.
-       *
-       * @see
-       * @ref GlossGhostedVector "vectors with ghost elements"
-       */
-      bool has_ghost_elements() const;
-
-      /**
-       * Return whether the vector contains only elements with value zero.
-       * This is a collective operation. This function is expensive, because
-       * potentially all elements have to be checked.
-       */
-      bool all_zero () const;
-
-      /**
-       * Return @p true if the vector has no negative entries, i.e. all
-       * entries are zero or positive. This function is used, for example, to
-       * check whether refinement indicators are really all positive (or
-       * zero).
-       *
-       * The function obviously only makes sense if the template argument of
-       * this class is a real type. If it is a complex type, then an exception
-       * is thrown.
-       */
-      bool is_non_negative () const;
-
-      /**
-       * Checks for equality of the two vectors.
-       */
-      template <typename Number2>
-      bool operator == (const Vector<Number2> &v) const;
-
-      /**
-       * Checks for inequality of the two vectors.
-       */
-      template <typename Number2>
-      bool operator != (const Vector<Number2> &v) const;
-
-      /**
-       * Perform the inner product of two vectors.
-       */
-      template <typename Number2>
-      Number operator * (const Vector<Number2> &V) const;
-
-      /**
-       * Computes the square of the l<sub>2</sub> norm of the vector (i.e.,
-       * the sum of the squares of all entries among all processors).
-       */
-      real_type norm_sqr () const;
-
-      /**
-       * Computes the mean value of all the entries in the vector.
-       */
-      Number mean_value () const;
-
-      /**
-       * Returns the l<sub>1</sub> norm of the vector (i.e., the sum of the
-       * absolute values of all entries among all processors).
-       */
-      real_type l1_norm () const;
-
-      /**
-       * Returns the l<sub>2</sub> norm of the vector (i.e., square root of
-       * the sum of the square of all entries among all processors).
-       */
-      real_type l2_norm () const;
-
-      /**
-       * Returns the l<sub>p</sub> norm with real @p p of the vector (i.e.,
-       * the pth root of sum of the pth power of all entries among all
-       * processors).
-       */
-      real_type lp_norm (const real_type p) const;
-
-      /**
-       * Returns the maximum norm of the vector (i.e., maximum absolute value
-       * among all entries among all processors).
-       */
-      real_type linfty_norm () const;
-
-      /**
-       * Performs a combined operation of a vector addition and a subsequent
-       * inner product, returning the value of the inner product. In other
-       * words, the result of this function is the same as if the user called
-       * @code
-       * this->add(a, V);
-       * return_value = *this * W;
-       * @endcode
-       *
-       * The reason this function exists is that this operation involves less
-       * memory transfer than calling the two functions separately. This
-       * method only needs to load three vectors, @p this, @p V, @p W, whereas
-       * calling separate methods means to load the calling vector @p this
-       * twice. Since most vector operations are memory transfer limited, this
-       * reduces the time by 25\% (or 50\% if @p W equals @p this).
-       */
-      Number add_and_dot (const Number          a,
-                          const Vector<Number> &V,
-                          const Vector<Number> &W);
-
-      /**
-       * Returns the global size of the vector, equal to the sum of the number
-       * of locally owned indices among all the processors.
-       */
-      size_type size () const;
-
-      /**
-       * Returns the local size of the vector, i.e., the number of indices
-       * owned locally.
-       */
-      size_type local_size() const;
-
-      /**
-       * Returns the half-open interval that specifies the locally owned range
-       * of the vector. Note that <code>local_size() == local_range().second -
-       * local_range().first</code>.
-       */
-      std::pair<size_type, size_type> local_range () const;
-
-      /**
-       * Returns true if the given global index is in the local range of this
-       * processor.
-       */
-      bool in_local_range (const size_type global_index) const;
-
-      /**
-       * Return an index set that describes which elements of this vector are
-       * owned by the current processor. Note that this index set does not
-       * include elements this vector may store locally as ghost elements but
-       * that are in fact owned by another processor. As a consequence, the
-       * index sets returned on different processors if this is a distributed
-       * vector will form disjoint sets that add up to the complete index set.
-       * Obviously, if a vector is created on only one processor, then the
-       * result would satisfy
-       * @code
-       *   vec.locally_owned_elements() == complete_index_set (vec.size())
-       * @endcode
-       */
-      IndexSet locally_owned_elements () const;
-
-      /**
-       * Returns the number of ghost elements present on the vector.
-       *
-       * This function is deprecated.
-       */
-      size_type n_ghost_entries () const DEAL_II_DEPRECATED;
-
-      /**
-       * Return an index set that describes which elements of this vector are
-       * not owned by the current processor but can be written into or read
-       * from locally (ghost elements).
-       *
-       * This function is deprecated.
-       */
-      const IndexSet &ghost_elements() const DEAL_II_DEPRECATED;
-
-      /**
-       * Returns whether the given global index is a ghost index on the
-       * present processor. Returns false for indices that are owned locally
-       * and for indices not present at all.
-       *
-       * This function is deprecated.
-       */
-      bool is_ghost_entry (const types::global_dof_index global_index) const DEAL_II_DEPRECATED;
-
-      /**
-       * Make the @p Vector class a bit like the <tt>vector<></tt> class of
-       * the C++ standard library by returning iterators to the start and end
-       * of the <i>locally owned</i> elements of this vector.
-       *
-       * It holds that end() - begin() == local_size().
-       */
-      iterator begin ();
-
-      /**
-       * Return constant iterator to the start of the locally owned elements
-       * of the vector.
-       */
-      const_iterator begin () const;
-
-      /**
-       * Return an iterator pointing to the element past the end of the array
-       * of locally owned entries.
-       */
-      iterator end ();
-
-      /**
-       * Return a constant iterator pointing to the element past the end of
-       * the array of the locally owned entries.
-       */
-      const_iterator end () const;
-      //@}
-
-
-      /**
-       * @name 2: Data-Access
-       */
-      //@{
-
-      /**
-       * Read access to the data in the position corresponding to @p
-       * global_index. The index must be either in the local range of the
-       * vector or be specified as a ghost index at construction.
-       *
-       * Performance: <tt>O(1)</tt> for locally owned elements that represent
-       * a contiguous range and <tt>O(log(n<sub>ranges</sub>))</tt> for ghost
-       * elements (quite fast, but slower than local_element()).
-       */
-      Number operator () (const size_type global_index) const;
-
-      /**
-       * Read and write access to the data in the position corresponding to @p
-       * global_index. The index must be either in the local range of the
-       * vector or be specified as a ghost index at construction.
-       *
-       * Performance: <tt>O(1)</tt> for locally owned elements that represent
-       * a contiguous range and <tt>O(log(n<sub>ranges</sub>))</tt> for ghost
-       * elements (quite fast, but slower than local_element()).
-       */
-      Number &operator () (const size_type global_index);
-
-      /**
-       * Read access to the data in the position corresponding to @p
-       * global_index. The index must be either in the local range of the
-       * vector or be specified as a ghost index at construction.
-       *
-       * This function does the same thing as operator().
-       */
-      Number operator [] (const size_type global_index) const;
-
-      /**
-       * Read and write access to the data in the position corresponding to @p
-       * global_index. The index must be either in the local range of the
-       * vector or be specified as a ghost index at construction.
-       *
-       * This function does the same thing as operator().
-       */
-      Number &operator [] (const size_type global_index);
-
-      /**
-       * A collective get operation: instead of getting individual elements of
-       * a vector, this function allows to get a whole set of elements at
-       * once. The indices of the elements to be read are stated in the first
-       * argument, the corresponding values are returned in the second.
-       */
-      template <typename OtherNumber>
-      void extract_subvector_to (const std::vector<size_type> &indices,
-                                 std::vector<OtherNumber> &values) const;
-
-      /**
-       * Just as the above, but with pointers. Useful in minimizing copying of
-       * data around.
-       */
-      template <typename ForwardIterator, typename OutputIterator>
-      void extract_subvector_to (ForwardIterator          indices_begin,
-                                 const ForwardIterator    indices_end,
-                                 OutputIterator           values_begin) const;
-
-      /**
-       * Read access to the data field specified by @p local_index. Locally
-       * owned indices can be accessed with indices
-       * <code>[0,local_size)</code>, and ghost indices with indices
-       * <code>[local_size,local_size+ n_ghost_entries]</code>.
-       *
-       * Performance: Direct array access (fast).
-       */
-      Number local_element (const size_type local_index) const;
-
-      /**
-       * Read and write access to the data field specified by @p local_index.
-       * Locally owned indices can be accessed with indices
-       * <code>[0,local_size)</code>, and ghost indices with indices
-       * <code>[local_size,local_size+n_ghosts]</code>.
-       *
-       * Performance: Direct array access (fast).
-       */
-      Number &local_element (const size_type local_index);
-      //@}
-
-
-      /**
-       * @name 3: Modification of vectors
-       */
-      //@{
-
-      /**
-       * Add the given vector to the present one.
-       */
-      Vector<Number> &operator += (const Vector<Number> &V);
-
-      /**
-       * Subtract the given vector from the present one.
-       */
-      Vector<Number> &operator -= (const Vector<Number> &V);
-
-      /**
-       * A collective add operation: This function adds a whole set of values
-       * stored in @p values to the vector components specified by @p indices.
-       */
-      template <typename OtherNumber>
-      void add (const std::vector<size_type>   &indices,
-                const std::vector<OtherNumber>  &values);
-
-      /**
-       * This is a second collective add operation. As a difference, this
-       * function takes a deal.II vector of values.
-       */
-      template <typename OtherNumber>
-      void add (const std::vector<size_type>        &indices,
-                const ::dealii::Vector<OtherNumber> &values);
-
-      /**
-       * Take an address where <tt>n_elements</tt> are stored contiguously and
-       * add them into the vector. Handles all cases which are not covered by
-       * the other two <tt>add()</tt> functions above.
-       */
-      template <typename OtherNumber>
-      void add (const size_type    n_elements,
-                const size_type   *indices,
-                const OtherNumber  *values);
-
-      /**
-       * Addition of @p s to all components. Note that @p s is a scalar and
-       * not a vector.
-       */
-      void add (const Number s);
-
-      /**
-       * Simple vector addition, equal to the <tt>operator +=</tt>.
-       *
-       * @deprecated Use the <tt>operator +=</tt> instead.
-       */
-      void add (const Vector<Number> &V) DEAL_II_DEPRECATED;
-
-      /**
-       * Simple addition of a multiple of a vector, i.e. <tt>*this +=
-       * a*V</tt>.
-       */
-      void add (const Number a, const Vector<Number> &V);
-
-      /**
-       * Multiple addition of scaled vectors, i.e. <tt>*this += a*V+b*W</tt>.
-       */
-      void add (const Number a, const Vector<Number> &V,
-                const Number b, const Vector<Number> &W);
-
-      /**
-       * Scaling and simple vector addition, i.e.  <tt>*this =
-       * s*(*this)+V</tt>.
-       */
-      void sadd (const Number          s,
-                 const Vector<Number> &V);
-
-      /**
-       * Scaling and simple addition, i.e.  <tt>*this = s*(*this)+a*V</tt>.
-       */
-      void sadd (const Number          s,
-                 const Number          a,
-                 const Vector<Number> &V);
-
-      /**
-       * Scaling and multiple addition.
-       *
-       * This function is deprecated.
-       */
-      void sadd (const Number          s,
-                 const Number          a,
-                 const Vector<Number> &V,
-                 const Number          b,
-                 const Vector<Number> &W) DEAL_II_DEPRECATED;
-
-      /**
-       * Scaling and multiple addition.  <tt>*this = s*(*this)+a*V + b*W +
-       * c*X</tt>.
-       *
-       * This function is deprecated.
-       */
-      void sadd (const Number          s,
-                 const Number          a,
-                 const Vector<Number> &V,
-                 const Number          b,
-                 const Vector<Number> &W,
-                 const Number          c,
-                 const Vector<Number> &X) DEAL_II_DEPRECATED;
-
-      /**
-       * Scale each element of the vector by a constant value.
-       */
-      Vector<Number> &operator *= (const Number factor);
-
-      /**
-       * Scale each element of the vector by the inverse of the given value.
-       */
-      Vector<Number> &operator /= (const Number factor);
-
-      /**
-       * Scale each element of this vector by the corresponding element in the
-       * argument. This function is mostly meant to simulate multiplication
-       * (and immediate re-assignment) by a diagonal scaling matrix.
-       */
-      void scale (const Vector<Number> &scaling_factors);
-
-      /**
-       * Scale each element of this vector by the corresponding element in the
-       * argument. This function is mostly meant to simulate multiplication
-       * (and immediate re-assignment) by a diagonal scaling matrix.
-       */
-      template <typename Number2>
-      void scale (const Vector<Number2> &scaling_factors);
-
-      /**
-       * Assignment <tt>*this = a*u</tt>.
-       */
-      void equ (const Number a, const Vector<Number> &u);
-
-      /**
-       * Assignment <tt>*this = a*u</tt>.
-       */
-      template <typename Number2>
-      void equ (const Number a, const Vector<Number2> &u);
-
-      /**
-       * Assignment <tt>*this = a*u + b*v</tt>.
-       *
-       * This function is deprecated.
-       */
-      void equ (const Number a, const Vector<Number> &u,
-                const Number b, const Vector<Number> &v) DEAL_II_DEPRECATED;
-
-      /**
-       * Assignment <tt>*this = a*u + b*v + b*w</tt>.
-       *
-       * This function is deprecated.
-       */
-      void equ (const Number a, const Vector<Number> &u,
-                const Number b, const Vector<Number> &v,
-                const Number c, const Vector<Number> &w) DEAL_II_DEPRECATED;
-
-      /**
-       * Compute the elementwise ratio of the two given vectors, that is let
-       * <tt>this[i] = a[i]/b[i]</tt>. This is useful for example if you want
-       * to compute the cellwise ratio of true to estimated error.
-       *
-       * This vector is appropriately scaled to hold the result.
-       *
-       * If any of the <tt>b[i]</tt> is zero, the result is undefined. No
-       * attempt is made to catch such situations.
-       */
-      void ratio (const Vector<Number> &a,
-                  const Vector<Number> &b) DEAL_II_DEPRECATED;
-      //@}
-
-
-      /**
-       * @name 4: Mixed stuff
-       */
-      //@{
-      /**
-       * Return a reference to the MPI communicator object in use with this
-       * vector.
-       */
-      const MPI_Comm &get_mpi_communicator () const;
-
-      /**
-       * Return the MPI partitioner that describes the parallel layout of the
-       * vector. This object can be used to initialize another vector with the
-       * respective reinit() call, for additional queries regarding the
-       * parallel communication, or the compatibility of partitioners.
-       */
-      std_cxx11::shared_ptr<const Utilities::MPI::Partitioner>
-      get_partitioner () const;
-
-      /**
-       * Checks whether the given partitioner is compatible with the
-       * partitioner used for this vector. Two partitioners are compatible if
-       * they have the same local size and the same ghost indices. They do not
-       * necessarily need to be the same data field of the shared pointer.
-       * This is a local operation only, i.e., if only some processors decide
-       * that the partitioning is not compatible, only these processors will
-       * return @p false, whereas the other processors will return @p true.
-       */
-      bool
-      partitioners_are_compatible (const Utilities::MPI::Partitioner &part) const;
-
-      /**
-       * Checks whether the given partitioner is compatible with the
-       * partitioner used for this vector. Two partitioners are compatible if
-       * they have the same local size and the same ghost indices. They do not
-       * necessarily need to be the same data field. As opposed to
-       * partitioners_are_compatible(), this method checks for compatibility
-       * among all processors and the method only returns @p true if the
-       * partitioner is the same on all processors.
-       *
-       * This method performs global communication, so make sure to use it
-       * only in a context where all processors call it the same number of
-       * times.
-       */
-      bool
-      partitioners_are_globally_compatible (const Utilities::MPI::Partitioner &part) const;
-
-      /**
-       * Prints the vector to the output stream @p out.
-       */
-      void print (std::ostream       &out,
-                  const unsigned int  precision  = 3,
-                  const bool          scientific = true,
-                  const bool          across     = true) const;
-
-      /**
-       * Returns the memory consumption of this class in bytes.
-       */
-      std::size_t memory_consumption () const;
-      //@}
-
-      /**
-       * Exception
-       */
-      DeclException3 (ExcNonMatchingElements,
-                      double, double, unsigned int,
-                      << "Called compress(VectorOperation::insert), but"
-                      << " the element received from a remote processor, value "
-                      << std::setprecision(16) << arg1
-                      << ", does not match with the value "
-                      << std::setprecision(16) << arg2
-                      << " on the owner processor " << arg3);
-
-      /**
-       * Exception
-       */
-      DeclException4 (ExcAccessToNonLocalElement,
-                      size_type, size_type, size_type, size_type,
-                      << "You tried to access element " << arg1
-                      << " of a distributed vector, but this element is not "
-                      << "stored on the current processor. Note: The range of "
-                      << "locally owned elements is " << arg2 << " to "
-                      << arg3 << ", and there are " << arg4 << " ghost elements "
-                      << "that this vector can access.");
-
-    private:
-      /**
-       * Local part of all_zero().
-       */
-      bool all_zero_local () const;
-
-      /**
-       * Local part of is_non_negative().
-       */
-      bool is_non_negative_local () const;
-
-      /**
-       * Local part of operator==.
-       */
-      template <typename Number2>
-      bool vectors_equal_local (const Vector<Number2> &v) const;
-
-      /**
-       * Local part of the inner product of two vectors.
-       */
-      template <typename Number2>
-      Number inner_product_local (const Vector<Number2> &V) const;
-
-      /**
-       * Local part of norm_sqr().
-       */
-      real_type norm_sqr_local () const;
-
-      /**
-       * Local part of mean_value().
-       */
-      Number mean_value_local () const;
-
-      /**
-       * Local part of l1_norm().
-       */
-      real_type l1_norm_local () const;
-
-      /**
-       * Local part of lp_norm().
-       */
-      real_type lp_norm_local (const real_type p) const;
-
-      /**
-       * Local part of linfty_norm().
-       */
-      real_type linfty_norm_local () const;
-
-      /**
-       * Local part of the addition followed by an inner product of two
-       * vectors.
-       */
-      Number add_and_dot_local (const Number          a,
-                                const Vector<Number> &V,
-                                const Vector<Number> &W);
-
-      /**
-       * Shared pointer to store the parallel partitioning information. This
-       * information can be shared between several vectors that have the same
-       * partitioning.
-       */
-      std_cxx11::shared_ptr<const Utilities::MPI::Partitioner> partitioner;
-
-      /**
-       * The size that is currently allocated in the val array.
-       */
-      size_type allocated_size;
-
-      /**
-       * Pointer to the array of local elements of this vector.
-       */
-      Number         *val;
-
-      /**
-       * Temporary storage that holds the data that is sent to this processor
-       * in @p compress() or sent from this processor in @p
-       * update_ghost_values.
-       */
-      mutable Number *import_data;
-
-      /**
-       * Stores whether the vector currently allows for reading ghost elements
-       * or not. Note that this is to ensure consistent ghost data and does
-       * not indicate whether the vector actually can store ghost elements. In
-       * particular, when assembling a vector we do not allow reading
-       * elements, only writing them.
-       */
-      mutable bool vector_is_ghosted;
-
-      /**
-       * Provide this class with all functionality of ::dealii::Vector by
-       * creating a VectorView object.
-       */
-      VectorView<Number> vector_view;
-
-#ifdef DEAL_II_WITH_MPI
-      /**
-       * A vector that collects all requests from @p compress() operations.
-       * This class uses persistent MPI communicators, i.e., the communication
-       * channels are stored during successive calls to a given function. This
-       * reduces the overhead involved with setting up the MPI machinery, but
-       * it does not remove the need for a receive operation to be posted
-       * before the data can actually be sent.
-       */
-      std::vector<MPI_Request>   compress_requests;
-
-      /**
-       * A vector that collects all requests from @p update_ghost_values()
-       * operations. This class uses persistent MPI communicators.
-       */
-      mutable std::vector<MPI_Request>   update_ghost_values_requests;
-#endif
-
-      /**
-       * A lock that makes sure that the @p compress and @p
-       * update_ghost_values functions give reasonable results also when used
-       * with several threads.
-       */
-      mutable Threads::Mutex mutex;
-
-      /**
-       * A helper function that clears the compress_requests and
-       * update_ghost_values_requests field. Used in reinit functions.
-       */
-      void clear_mpi_requests ();
-
-      /**
-       * A helper function that is used to resize the val array.
-       */
-      void resize_val (const size_type new_allocated_size);
-
-      /*
-       * Make all other vector types friends.
-       */
-      template <typename Number2> friend class Vector;
-
-      /**
-       * Make BlockVector type friends.
-       */
-      template <typename Number2> friend class BlockVector;
-    };
-
-    /*@}*/
-
-
-    /*----------------------- Inline functions ----------------------------------*/
-
-#ifndef DOXYGEN
-
-    template <typename Number>
-    inline
-    Vector<Number>::Vector ()
-      :
-      partitioner (new Utilities::MPI::Partitioner()),
-      allocated_size (0),
-      val (0),
-      import_data (0),
-      vector_is_ghosted (false),
-      vector_view (0, static_cast<Number *>(0))
-    {}
-
-
-
-    template <typename Number>
-    inline
-    Vector<Number>::Vector (const Vector<Number> &v)
-      :
-      Subscriptor(),
-      allocated_size (0),
-      val (0),
-      import_data (0),
-      vector_is_ghosted (false),
-      vector_view (0, static_cast<Number *>(0))
-    {
-      reinit (v, true);
-      vector_view = v.vector_view;
-      zero_out_ghosts();
-    }
-
-
-
-    template <typename Number>
-    inline
-    Vector<Number>::Vector (const IndexSet &local_range,
-                            const IndexSet &ghost_indices,
-                            const MPI_Comm  communicator)
-      :
-      allocated_size (0),
-      val (0),
-      import_data (0),
-      vector_is_ghosted (false),
-      vector_view (0, static_cast<Number *>(0))
-    {
-      reinit (local_range, ghost_indices, communicator);
-    }
-
-
-
-    template <typename Number>
-    inline
-    Vector<Number>::Vector (const IndexSet &local_range,
-                            const MPI_Comm  communicator)
-      :
-      allocated_size (0),
-      val (0),
-      import_data (0),
-      vector_is_ghosted (false),
-      vector_view (0, static_cast<Number *>(0))
-    {
-      reinit (local_range, communicator);
-    }
-
-
-
-    template <typename Number>
-    inline
-    Vector<Number>::Vector (const size_type size)
-      :
-      allocated_size (0),
-      val (0),
-      import_data (0),
-      vector_is_ghosted (false),
-      vector_view (0, static_cast<Number *>(0))
-    {
-      reinit (size, false);
-    }
-
-
-
-    template <typename Number>
-    inline
-    Vector<Number>::
-    Vector (const std_cxx11::shared_ptr<const Utilities::MPI::Partitioner> &partitioner)
-      :
-      allocated_size (0),
-      val (0),
-      import_data (0),
-      vector_is_ghosted (false),
-      vector_view (0, static_cast<Number *>(0))
-    {
-      reinit (partitioner);
-    }
-
-
-
-    template <typename Number>
-    inline
-    Vector<Number>::~Vector ()
-    {
-      resize_val(0);
-
-      if (import_data != 0)
-        delete[] import_data;
-      import_data = 0;
-
-      clear_mpi_requests();
-    }
-
-
-
-    template <typename Number>
-    inline
-    Vector<Number> &
-    Vector<Number>::operator = (const Vector<Number> &c)
-    {
-#ifdef _MSC_VER
-      return this->operator=<Number>(c);
-#else
-      return this->template operator=<Number>(c);
-#endif
-    }
-
-
-
-    template <typename Number>
-    template <typename Number2>
-    inline
-    Vector<Number> &
-    Vector<Number>::operator = (const Vector<Number2> &c)
-    {
-      Assert (c.partitioner.get() != 0, ExcNotInitialized());
-
-      // we update ghost values whenever one of the input or output vector
-      // already held ghost values or when we import data from a vector with
-      // the same local range but different ghost layout
-      bool must_update_ghost_values = c.vector_is_ghosted;
-
-      // check whether the two vectors use the same parallel partitioner. if
-      // not, check if all local ranges are the same (that way, we can
-      // exchange data between different parallel layouts). One variant which
-      // is included here and necessary for compatibility with the other
-      // distributed vector classes (Trilinos, PETSc) is the case when vector
-      // c does not have any ghosts (constructed without ghost elements given)
-      // but the current vector does: In that case, we need to exchange data
-      // also when none of the two vector had updated its ghost values before.
-      if (partitioner.get() == 0)
-        reinit (c, true);
-      else if (partitioner.get() != c.partitioner.get())
-        {
-          // local ranges are also the same if both partitioners are empty
-          // (even if they happen to define the empty range as [0,0) or [c,c)
-          // for some c!=0 in a different way).
-          int local_ranges_are_identical =
-            (local_range() == c.local_range() ||
-             (local_range().second == local_range().first &&
-              c.local_range().second == c.local_range().first));
-          if ((c.partitioner->n_mpi_processes() > 1 &&
-               Utilities::MPI::min(local_ranges_are_identical,
-                                   c.partitioner->get_communicator()) == 0)
-              ||
-              !local_ranges_are_identical)
-            reinit (c, true);
-          else
-            must_update_ghost_values |= vector_is_ghosted;
-
-          must_update_ghost_values |=
-            (c.partitioner->ghost_indices_initialized() == false &&
-             partitioner->ghost_indices_initialized() == true);
-        }
-      else
-        must_update_ghost_values |= vector_is_ghosted;
-
-      // Need to explicitly downcast to dealii::Vector to make templated
-      // operator= available.
-      AssertDimension(vector_view.size(), c.vector_view.size());
-      static_cast<dealii::Vector<Number> &>(vector_view) = c.vector_view;
-
-      if (must_update_ghost_values)
-        update_ghost_values();
-      else
-        zero_out_ghosts();
-      return *this;
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::compress (::dealii::VectorOperation::values operation)
-    {
-      compress_start (0, operation);
-      compress_finish(operation);
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::update_ghost_values () const
-    {
-      update_ghost_values_start ();
-      update_ghost_values_finish ();
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::zero_out_ghosts ()
-    {
-      std::fill_n (&val[partitioner->local_size()],
-                   partitioner->n_ghost_indices(),
-                   Number());
-      vector_is_ghosted = false;
-    }
-
-
-
-    template <typename Number>
-    inline
-    bool
-    Vector<Number>::has_ghost_elements () const
-    {
-      return vector_is_ghosted;
-    }
-
-
-
-    template <typename Number>
-    inline
-    bool
-    Vector<Number>::all_zero_local () const
-    {
-      return partitioner->local_size()>0 ? vector_view.all_zero () : true;
-    }
-
-
-
-    template <typename Number>
-    inline
-    bool
-    Vector<Number>::all_zero () const
-    {
-      // use int instead of bool. in order to make global reduction operations
-      // work also when MPI_Init was not called, only call MPI_Allreduce
-      // commands when there is more than one processor (note that reinit()
-      // functions handle this case correctly through the job_supports_mpi()
-      // query). this is the same in all the functions below
-      int local_result = -static_cast<int>(all_zero_local());
-      if (partitioner->n_mpi_processes() > 1)
-        return -Utilities::MPI::max(local_result,
-                                    partitioner->get_communicator());
-      else
-        return -local_result;
-    }
-
-
-
-    template <typename Number>
-    inline
-    bool
-    Vector<Number>::is_non_negative_local () const
-    {
-      return partitioner->local_size()>0 ? vector_view.is_non_negative () : true;
-    }
-
-
-
-    template <typename Number>
-    inline
-    bool
-    Vector<Number>::is_non_negative () const
-    {
-      int local_result = -static_cast<int>(is_non_negative_local());
-      if (partitioner->n_mpi_processes() > 1)
-        return -Utilities::MPI::max(local_result,
-                                    partitioner->get_communicator());
-      else
-        return -local_result;
-    }
-
-
-
-    template <typename Number>
-    template <typename Number2>
-    inline
-    bool
-    Vector<Number>::vectors_equal_local (const Vector<Number2> &v) const
-    {
-      return partitioner->local_size()>0 ?
-             vector_view.template operator == <Number2>(v.vector_view)
-             : true;
-    }
-
-
-
-    template <typename Number>
-    template <typename Number2>
-    inline
-    bool
-    Vector<Number>::operator == (const Vector<Number2> &v) const
-    {
-      // MPI does not support bools, so use unsigned int instead. Two vectors
-      // are equal if the check for non-equal fails on all processors
-      unsigned int local_result = static_cast<int>(!vectors_equal_local(v));
-      unsigned int result =
-        partitioner->n_mpi_processes() > 1
-        ?
-        Utilities::MPI::max(local_result, partitioner->get_communicator())
-        :
-        local_result;
-      return result==0;
-    }
-
-
-
-    template <typename Number>
-    template <typename Number2>
-    inline
-    bool
-    Vector<Number>::operator != (const Vector<Number2> &v) const
-    {
-      return !(operator == (v));
-    }
-
-
-
-    template <typename Number>
-    template <typename Number2>
-    inline
-    Number
-    Vector<Number>::inner_product_local(const Vector<Number2> &V) const
-    {
-      // on some processors, the size might be zero, which is not allowed by
-      // the dealii::Vector class. Therefore, insert a check here
-      return (partitioner->local_size()>0 ?
-              vector_view.operator* (V.vector_view)
-              : Number());
-    }
-
-
-
-    template <typename Number>
-    template <typename Number2>
-    inline
-    Number
-    Vector<Number>::operator * (const Vector<Number2> &V) const
-    {
-      Number local_result = inner_product_local(V);
-      if (partitioner->n_mpi_processes() > 1)
-        return Utilities::MPI::sum (local_result,
-                                    partitioner->get_communicator());
-      else
-        return local_result;
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::real_type
-    Vector<Number>::norm_sqr_local () const
-    {
-      return partitioner->local_size()>0 ? vector_view.norm_sqr() : real_type();
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::real_type
-    Vector<Number>::norm_sqr () const
-    {
-      real_type local_result = norm_sqr_local();
-      if (partitioner->n_mpi_processes() > 1)
-        return Utilities::MPI::sum(local_result,
-                                   partitioner->get_communicator());
-      else
-        return local_result;
-    }
-
-
-
-    template <typename Number>
-    inline
-    Number
-    Vector<Number>::mean_value_local () const
-    {
-      Assert (partitioner->size()!=0, ExcEmptyObject());
-      return (partitioner->local_size() ?
-              vector_view.mean_value()
-              : Number());
-    }
-
-
-
-    template <typename Number>
-    inline
-    Number
-    Vector<Number>::mean_value () const
-    {
-      Number local_result = mean_value_local();
-      if (partitioner->n_mpi_processes() > 1)
-        return Utilities::MPI::sum (local_result *
-                                    (real_type)partitioner->local_size(),
-                                    partitioner->get_communicator())
-               /(real_type)partitioner->size();
-      else
-        return local_result;
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::real_type
-    Vector<Number>::l1_norm_local () const
-    {
-      return partitioner->local_size() ? vector_view.l1_norm() : real_type();
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::real_type
-    Vector<Number>::l1_norm () const
-    {
-      real_type local_result = l1_norm_local();
-      if (partitioner->n_mpi_processes() > 1)
-        return Utilities::MPI::sum(local_result,
-                                   partitioner->get_communicator());
-      else
-        return local_result;
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::real_type
-    Vector<Number>::l2_norm () const
-    {
-      return std::sqrt(norm_sqr());
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::real_type
-    Vector<Number>::lp_norm_local (const real_type p) const
-    {
-      return partitioner->local_size() ? vector_view.lp_norm(p) : real_type();
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::real_type
-    Vector<Number>::lp_norm (const real_type p) const
-    {
-      const real_type local_result = lp_norm_local(p);
-      if (partitioner->n_mpi_processes() > 1)
-        return std::pow (Utilities::MPI::sum(std::pow(local_result,p),
-                                             partitioner->get_communicator()),
-                         static_cast<real_type>(1.0/p));
-      else
-        return local_result;
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::real_type
-    Vector<Number>::linfty_norm_local () const
-    {
-      return partitioner->local_size() ? vector_view.linfty_norm() : real_type();
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::real_type
-    Vector<Number>::linfty_norm () const
-    {
-      const real_type local_result = linfty_norm_local();
-      if (partitioner->n_mpi_processes() > 1)
-        return Utilities::MPI::max (local_result,
-                                    partitioner->get_communicator());
-      else
-        return local_result;
-    }
-
-
-
-    template <typename Number>
-    inline
-    Number
-    Vector<Number>::add_and_dot_local(const Number          a,
-                                      const Vector<Number> &V,
-                                      const Vector<Number> &W)
-    {
-      // on some processors, the size might be zero, which is not allowed by
-      // the dealii::Vector class. Therefore, insert a check here
-      return (partitioner->local_size()>0 ?
-              vector_view.add_and_dot(a, V.vector_view, W.vector_view)
-              : Number());
-    }
-
-
-
-    template <typename Number>
-    inline
-    Number
-    Vector<Number>::add_and_dot (const Number          a,
-                                 const Vector<Number> &V,
-                                 const Vector<Number> &W)
-    {
-      Number local_result = add_and_dot_local(a, V, W);
-      if (partitioner->n_mpi_processes() > 1)
-        return Utilities::MPI::sum (local_result,
-                                    partitioner->get_communicator());
-      else
-        return local_result;
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::size_type
-    Vector<Number>::size () const
-    {
-      return partitioner->size();
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::size_type
-    Vector<Number>::local_size () const
-    {
-      return partitioner->local_size();
-    }
-
-
-
-    template <typename Number>
-    inline
-    std::pair<typename Vector<Number>::size_type,
-        typename Vector<Number>::size_type>
-        Vector<Number>::local_range () const
-    {
-      return partitioner->local_range();
-    }
-
-
-
-    template <typename Number>
-    inline
-    bool
-    Vector<Number>::in_local_range
-    (const size_type global_index) const
-    {
-      return partitioner->in_local_range (global_index);
-    }
-
-
-
-    template <typename Number>
-    inline
-    IndexSet
-    Vector<Number>::locally_owned_elements() const
-    {
-      IndexSet is (size());
-
-      is.add_range (local_range().first, local_range().second);
-
-      return is;
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::size_type
-    Vector<Number>::n_ghost_entries () const
-    {
-      return partitioner->n_ghost_indices();
-    }
-
-
-
-    template <typename Number>
-    inline
-    const IndexSet &
-    Vector<Number>::ghost_elements() const
-    {
-      return partitioner->ghost_indices();
-    }
-
-
-
-    template <typename Number>
-    inline
-    bool
-    Vector<Number>::is_ghost_entry (const size_type global_index) const
-    {
-      return partitioner->is_ghost_entry (global_index);
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::iterator
-    Vector<Number>::begin ()
-    {
-      return vector_view.begin();
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::const_iterator
-    Vector<Number>::begin () const
-    {
-      return vector_view.begin();
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::iterator
-    Vector<Number>::end ()
-    {
-      return vector_view.end();
-    }
-
-
-
-    template <typename Number>
-    inline
-    typename Vector<Number>::const_iterator
-    Vector<Number>::end () const
-    {
-      return vector_view.end();
-    }
-
-
-
-    template <typename Number>
-    inline
-    Number
-    Vector<Number>::operator() (const size_type global_index) const
-    {
-      Assert (in_local_range (global_index) ||
-              partitioner->ghost_indices().is_element(global_index),
-              ExcAccessToNonLocalElement(global_index, local_range().first,
-                                         local_range().second,
-                                         partitioner->ghost_indices().n_elements()));
-      // do not allow reading a vector which is not in ghost mode
-      Assert (in_local_range (global_index) || vector_is_ghosted == true,
-              ExcMessage("You tried to read a ghost element of this vector, "
-                         "but it has not imported its ghost values."));
-      return val[partitioner->global_to_local(global_index)];
-    }
-
-
-
-    template <typename Number>
-    inline
-    Number &
-    Vector<Number>::operator() (const size_type global_index)
-    {
-      Assert (in_local_range (global_index) ||
-              partitioner->ghost_indices().is_element(global_index),
-              ExcAccessToNonLocalElement(global_index, local_range().first,
-                                         local_range().second,
-                                         partitioner->ghost_indices().n_elements()));
-      // we would like to prevent reading ghosts from a vector that does not
-      // have them imported, but this is not possible because we might be in a
-      // part of the code where the vector has enabled ghosts but is non-const
-      // (then, the compiler picks this method according to the C++ rule book
-      // even if a human would pick the const method when this subsequent use
-      // is just a read)
-      return val[partitioner->global_to_local (global_index)];
-    }
-
-
-
-    template <typename Number>
-    inline
-    Number
-    Vector<Number>::operator[] (const size_type global_index) const
-    {
-      return operator()(global_index);
-    }
-
-
-
-    template <typename Number>
-    inline
-    Number &
-    Vector<Number>::operator[] (const size_type global_index)
-    {
-      return operator()(global_index);
-    }
-
-
-
-    template <typename Number>
-    template <typename OtherNumber>
-    inline
-    void Vector<Number>::extract_subvector_to (const std::vector<size_type> &indices,
-                                               std::vector<OtherNumber> &values) const
-    {
-      for (size_type i = 0; i < indices.size(); ++i)
-        values[i] = operator()(indices[i]);
-    }
-
-
-
-    template <typename Number>
-    template <typename ForwardIterator, typename OutputIterator>
-    inline
-    void Vector<Number>::extract_subvector_to (ForwardIterator          indices_begin,
-                                               const ForwardIterator    indices_end,
-                                               OutputIterator           values_begin) const
-    {
-      while (indices_begin != indices_end)
-        {
-          *values_begin = operator()(*indices_begin);
-          indices_begin++;
-          values_begin++;
-        }
-    }
-
-
-
-    template <typename Number>
-    inline
-    Number
-    Vector<Number>::local_element (const size_type local_index) const
-    {
-      AssertIndexRange (local_index,
-                        partitioner->local_size()+
-                        partitioner->n_ghost_indices());
-      // do not allow reading a vector which is not in ghost mode
-      Assert (local_index < local_size() || vector_is_ghosted == true,
-              ExcMessage("You tried to read a ghost element of this vector, "
-                         "but it has not imported its ghost values."));
-      return val[local_index];
-    }
-
-
-
-    template <typename Number>
-    inline
-    Number &
-    Vector<Number>::local_element (const size_type local_index)
-    {
-      AssertIndexRange (local_index,
-                        partitioner->local_size()+
-                        partitioner->n_ghost_indices());
-      return val[local_index];
-    }
-
-
-
-    template <typename Number>
-    inline
-    Vector<Number> &
-    Vector<Number>::operator = (const Number s)
-    {
-      // if we call Vector::operator=0, we want to zero out all the entries
-      // plus ghosts.
-      if (partitioner->local_size() > 0)
-        vector_view.dealii::template Vector<Number>::operator= (s);
-      if (s==Number())
-        zero_out_ghosts();
-
-      return *this;
-    }
-
-
-
-    template <typename Number>
-    inline
-    Vector<Number> &
-    Vector<Number>::operator += (const Vector<Number> &v)
-    {
-      AssertDimension (local_size(), v.local_size());
-
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size()>0)
-        vector_view += v.vector_view;
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-
-      return *this;
-    }
-
-
-
-    template <typename Number>
-    inline
-    Vector<Number> &
-    Vector<Number>::operator -= (const Vector<Number> &v)
-    {
-      AssertDimension (local_size(), v.local_size());
-
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size()>0)
-        vector_view -= v.vector_view;
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-
-      return *this;
-    }
-
-
-
-    template <typename Number>
-    template <typename OtherNumber>
-    inline
-    void
-    Vector<Number>::add (const std::vector<size_type> &indices,
-                         const std::vector<OtherNumber>  &values)
-    {
-      AssertDimension (indices.size(), values.size());
-      add (indices.size(), &indices[0], &values[0]);
-    }
-
-
-
-    template <typename Number>
-    template <typename OtherNumber>
-    inline
-    void
-    Vector<Number>::add (const std::vector<size_type>    &indices,
-                         const ::dealii::Vector<OtherNumber> &values)
-    {
-      AssertDimension (indices.size(), values.size());
-      add (indices.size(), &indices[0], values.begin());
-    }
-
-
-
-    template <typename Number>
-    template <typename OtherNumber>
-    inline
-    void
-    Vector<Number>::add (const size_type    n_indices,
-                         const size_type   *indices,
-                         const OtherNumber *values)
-    {
-      for (size_type i=0; i<n_indices; ++i)
-        {
-          Assert (numbers::is_finite(values[i]),
-                  ExcMessage("The given value is not finite but either infinite or Not A Number (NaN)"));
-          this->operator()(indices[i]) += values[i];
-        }
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::add (const Number a)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.add (a);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::add (const Vector<Number> &v)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.add (v.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::add (const Number a,
-                         const Vector<Number> &v)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.add (a, v.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::add (const Number a,
-                         const Vector<Number> &v,
-                         const Number b,
-                         const Vector<Number> &w)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.add (a, v.vector_view, b, w.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::sadd (const Number x,
-                          const Vector<Number> &v)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.sadd (x, v.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::sadd (const Number x,
-                          const Number a,
-                          const Vector<Number> &v)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.sadd (x, a, v.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::sadd (const Number x,
-                          const Number a,
-                          const Vector<Number> &v,
-                          const Number b,
-                          const Vector<Number> &w)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.sadd (x, a, v.vector_view, b, w.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::sadd (const Number s,
-                          const Number a,
-                          const Vector<Number> &v,
-                          const Number b,
-                          const Vector<Number> &w,
-                          const Number c,
-                          const Vector<Number> &x)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.sadd (s, a, v.vector_view, b, w.vector_view,
-                          c, x.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    inline
-    Vector<Number> &
-    Vector<Number>::operator *= (const Number factor)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view *= factor;
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-
-      return *this;
-    }
-
-
-
-    template <typename Number>
-    inline
-    Vector<Number> &
-    Vector<Number>::operator /= (const Number factor)
-    {
-      operator *= (1./factor);
-      return *this;
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::scale (const Vector<Number> &scaling_factors)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.scale (scaling_factors.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    template <typename Number2>
-    inline
-    void
-    Vector<Number>::scale (const Vector<Number2> &scaling_factors)
-    {
-      if (local_size())
-        vector_view.template scale<Number2> (scaling_factors.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::equ (const Number a,
-                         const Vector<Number> &v)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.equ (a, v.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    template <typename Number2>
-    inline
-    void
-    Vector<Number>::equ (const Number a,
-                         const Vector<Number2> &v)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.equ (a, v.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::equ (const Number a,
-                         const Vector<Number> &v,
-                         const Number b,
-                         const Vector<Number> &w)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.equ (a, v.vector_view, b, w.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::equ (const Number a,
-                         const Vector<Number> &v,
-                         const Number b,
-                         const Vector<Number> &w,
-                         const Number c,
-                         const Vector<Number> &x)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.equ (a, v.vector_view, b, w.vector_view,
-                         c, x.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    inline
-    void
-    Vector<Number>::ratio (const Vector<Number> &a,
-                           const Vector<Number> &b)
-    {
-      // dealii::Vector does not allow empty fields but this might happen on
-      // some processors for parallel implementation
-      if (local_size())
-        vector_view.ratio (a.vector_view, b.vector_view);
-
-      if (vector_is_ghosted)
-        update_ghost_values();
-    }
-
-
-
-    template <typename Number>
-    inline
-    const MPI_Comm &
-    Vector<Number>::get_mpi_communicator() const
-    {
-      return partitioner->get_communicator();
-    }
-
-
-
-    template <typename Number>
-    inline
-    std_cxx11::shared_ptr<const Utilities::MPI::Partitioner>
-    Vector<Number>::get_partitioner () const
-    {
-      return partitioner;
-    }
-
-#endif  // ifndef DOXYGEN
-
-  } // end of namespace distributed
-
-} // end of namespace parallel
-
-
-
-/**
- * Global function @p swap which overloads the default implementation of the
- * C++ standard library which uses a temporary object. The function simply
- * exchanges the data of the two vectors.
- *
- * @relates Vector
- * @author Katharina Kormann, Martin Kronbichler, 2011
- */
-template <typename Number>
-inline
-void swap (parallel::distributed::Vector<Number> &u,
-           parallel::distributed::Vector<Number> &v)
-{
-  u.swap (v);
+    using LinearAlgebra::parallel::Vector;
+  }
 }
 
-
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
index 5cf413b32d1edf01981abb487635e667c1e6408b..b633050b42a9a5b24417c541492dbd4e050693a6 100644 (file)
@@ -33,9 +33,9 @@ DEAL_II_NAMESPACE_OPEN
 
 template <typename number> class Vector;
 template <typename number> class SparseMatrix;
-namespace parallel
+namespace LinearAlgebra
 {
-  namespace distributed
+  namespace parallel
   {
     template <typename number> class Vector;
   }
@@ -1715,14 +1715,14 @@ namespace internal
     template <typename Number>
     inline
     void
-    vector_updates (const parallel::distributed::Vector<Number> &src,
-                    const parallel::distributed::Vector<Number> &matrix_diagonal_inverse,
+    vector_updates (const LinearAlgebra::parallel::Vector<Number> &src,
+                    const LinearAlgebra::parallel::Vector<Number> &matrix_diagonal_inverse,
                     const bool    start_zero,
                     const double  factor1,
                     const double  factor2,
-                    parallel::distributed::Vector<Number> &update1,
-                    parallel::distributed::Vector<Number> &update2,
-                    parallel::distributed::Vector<Number> &dst)
+                    LinearAlgebra::parallel::Vector<Number> &update1,
+                    LinearAlgebra::parallel::Vector<Number> &update2,
+                    LinearAlgebra::parallel::Vector<Number> &dst)
     {
       VectorUpdater<Number> upd(src.begin(), matrix_diagonal_inverse.begin(),
                                 start_zero, factor1, factor2,
index 00c1d7d31d14c18fbf0b8a7bbd7861ccdfd0826c..d42be66c2f18759110c9c79372676c2c9648eff6 100644 (file)
@@ -529,7 +529,7 @@ namespace LinearAlgebra
      * For parallel loops with TBB, this member variable stores the affinity
      * information of loops.
      */
-    mutable std_cxx11::shared_ptr<parallel::internal::TBBPartitioner> thread_loop_partitioner;
+    mutable std_cxx11::shared_ptr<::dealii::parallel::internal::TBBPartitioner> thread_loop_partitioner;
 
     /**
      * Make all other ReadWriteVector types friends.
index 7fadf02cda4be1e332e2944160e8b831917fdaa9..593256c0ebb83f9eb85cc327991605dc10ae1506 100644 (file)
@@ -902,6 +902,11 @@ namespace internal
         // order to obtain known loop bounds for most of the work.
         size_type index = first;
         ResultType outer_results [vector_accumulation_recursion_threshold];
+
+        // set the zeroth element to zero to correctly handle the case where
+        // vec_size == 0
+        outer_results[0] = ResultType();
+
         size_type n_chunks = vec_size / 32;
         const size_type remainder = vec_size % 32;
         Assert (remainder == 0 || n_chunks < vector_accumulation_recursion_threshold,
@@ -1214,7 +1219,7 @@ namespace internal
         partitioner->release_one_partitioner(tbb_partitioner);
         result = generic_functor.do_sum();
       }
-    else if (vec_size > 0)
+    else
       accumulate_recursive(op,0,vec_size,result);
 #else
     accumulate_recursive(op,0,vec_size,result);
index 5decb1574e73764f0019295e0085d39322601d8b..de05d3e1e4f32e5d66ae72da8c655867063eb684 100644 (file)
 
 DEAL_II_NAMESPACE_OPEN
 
-class CommunicationPatternBase;
 class IndexSet;
-template <typename Number> class ReadWriteVector;
+namespace LinearAlgebra
+{
+  class CommunicationPatternBase;
+  template <typename Number> class ReadWriteVector;
+}
 
 namespace LinearAlgebra
 {
index b8dcb9fc3aeb8a608cf145261f075e3906aafa56..0903bc2f377a324502f8d3df95182e7531bb644c 100644 (file)
@@ -34,9 +34,9 @@ DEAL_II_NAMESPACE_OPEN
 
 
 // forward declarations
-namespace parallel
+namespace LinearAlgebra
 {
-  namespace distributed
+  namespace parallel
   {
     template <typename> class Vector;
   }
@@ -2117,7 +2117,7 @@ namespace internal
   template <typename Number>
   inline
   Number &
-  vector_access (parallel::distributed::Vector<Number> &vec,
+  vector_access (LinearAlgebra::parallel::Vector<Number> &vec,
                  const unsigned int                     entry)
   {
     return vec.local_element(entry);
@@ -2131,7 +2131,7 @@ namespace internal
   template <typename Number>
   inline
   Number
-  vector_access (const parallel::distributed::Vector<Number> &vec,
+  vector_access (const LinearAlgebra::parallel::Vector<Number> &vec,
                  const unsigned int                           entry)
   {
     return vec.local_element(entry);
@@ -2155,7 +2155,7 @@ namespace internal
 
   template <typename Number>
   inline
-  void check_vector_compatibility (const parallel::distributed::Vector<Number>  &vec,
+  void check_vector_compatibility (const LinearAlgebra::parallel::Vector<Number>  &vec,
                                    const internal::MatrixFreeFunctions::DoFInfo &dof_info)
   {
     Assert (vec.partitioners_are_compatible(*dof_info.vector_partitioner),
index 3e44d8b2a31132acb9c501336fdf75987f4a732b..abeae873e4aff05462c5dac065e6512d9f614a07 100644 (file)
@@ -89,6 +89,19 @@ namespace Utilities
 
 
 
+    void
+    Partitioner::reinit(const IndexSet &vector_space_vector_index_set,
+                        const IndexSet &read_write_vector_index_set,
+                        const MPI_Comm &communicator_in)
+    {
+      have_ghost_indices = false;
+      communicator = communicator_in;
+      set_owned_indices (vector_space_vector_index_set);
+      set_ghost_indices (read_write_vector_index_set);
+    }
+
+
+
     void
     Partitioner::set_owned_indices (const IndexSet &locally_owned_indices)
     {
index 9dc145c2c5943827853d36622409fa4998343262..9ead26e3fef7a520500424a64f8a896f2b4f37b6 100644 (file)
@@ -28,9 +28,9 @@ SET(_src
   full_matrix.cc
   lapack_full_matrix.cc
   la_vector.cc
+  la_parallel_vector.cc
   matrix_lib.cc
   matrix_out.cc
-  parallel_vector.cc
   precondition_block.cc
   precondition_block_ez.cc
   relaxation_block.cc
@@ -62,7 +62,7 @@ SET(_inst
   full_matrix.inst.in
   lapack_full_matrix.inst.in
   la_vector.inst.in
-  parallel_vector.inst.in
+  la_parallel_vector.inst.in
   precondition_block.inst.in
   relaxation_block.inst.in
   read_write_vector.inst.in
diff --git a/source/lac/la_parallel_vector.cc b/source/lac/la_parallel_vector.cc
new file mode 100644 (file)
index 0000000..f62ac75
--- /dev/null
@@ -0,0 +1,50 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/lac/la_parallel_vector.h>
+#include <deal.II/lac/la_parallel_vector.templates.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+#include "la_parallel_vector.inst"
+
+// do a few functions that currently don't fit the scheme because they have
+// two template arguments that need to be different (the case of same
+// arguments is covered by the default copy constructor and copy operator that
+// is declared separately)
+
+namespace LinearAlgebra
+{
+  namespace parallel
+  {
+#define TEMPL_COPY_CONSTRUCTOR(S1,S2)                                   \
+  template Vector<S1>& Vector<S1>::operator=<S2> (const Vector<S2> &)
+
+#ifndef DEAL_II_EXPLICIT_CONSTRUCTOR_BUG
+    TEMPL_COPY_CONSTRUCTOR(double,float);
+    TEMPL_COPY_CONSTRUCTOR(float,double);
+
+
+    TEMPL_COPY_CONSTRUCTOR(std::complex<double>,std::complex<float>);
+    TEMPL_COPY_CONSTRUCTOR(std::complex<float>,std::complex<double>);
+
+#endif
+
+#undef TEMPL_COPY_CONSTRUCTOR
+  }
+}
+
+
+DEAL_II_NAMESPACE_CLOSE
similarity index 90%
rename from source/lac/parallel_vector.inst.in
rename to source/lac/la_parallel_vector.inst.in
index 1291f9eeeb5d8ba1dd1dad3dfacda8503cf95762..99b7760bd4fb2545bdce676db619450dce5354d9 100644 (file)
@@ -17,9 +17,9 @@
 
 for (SCALAR : REAL_SCALARS)
 {
-  namespace parallel
+  namespace LinearAlgebra
   \{
-    namespace distributed
+    namespace parallel
     \{
       template class Vector<SCALAR>;
     \}
@@ -28,9 +28,9 @@ for (SCALAR : REAL_SCALARS)
 
 for (S1, S2 : REAL_SCALARS)
 {
-  namespace parallel
+  namespace LinearAlgebra
   \{
-    namespace distributed
+    namespace parallel
     \{
       template void Vector<S1>::reinit<S2> (const Vector<S2>&,
                                             const bool);
diff --git a/source/lac/parallel_vector.cc b/source/lac/parallel_vector.cc
deleted file mode 100644 (file)
index fd01c6d..0000000
+++ /dev/null
@@ -1,24 +0,0 @@
-// ---------------------------------------------------------------------
-//
-// Copyright (C) 1999 - 2014 by the deal.II authors
-//
-// This file is part of the deal.II library.
-//
-// The deal.II library is free software; you can use it, redistribute
-// it, and/or modify it under the terms of the GNU Lesser General
-// Public License as published by the Free Software Foundation; either
-// version 2.1 of the License, or (at your option) any later version.
-// The full text of the license can be found in the file LICENSE at
-// the top level of the deal.II distribution.
-//
-// ---------------------------------------------------------------------
-
-#include <deal.II/lac/parallel_vector.h>
-#include <deal.II/lac/parallel_vector.templates.h>
-
-DEAL_II_NAMESPACE_OPEN
-
-#include "parallel_vector.inst"
-
-
-DEAL_II_NAMESPACE_CLOSE
index 2460d8226b0d72f58c4f028dd03e59524f2a1637..0ad8c1e430eb2566b742716142900a25d38a3f2f 100644 (file)
@@ -105,7 +105,7 @@ void test ()
   }
   // check inner product
   {
-    const double norm_sqr = v.norm_sqr();
+    const double norm_sqr = v.l2_norm() * v.l2_norm();
     AssertThrow (std::fabs(v * v - norm_sqr) < 1e-15,
                  ExcInternalError());
     parallel::distributed::Vector<double> v2;
@@ -120,48 +120,6 @@ void test ()
       deallog << "Inner product: " << inner_prod << std::endl;
   }
 
-  // check operator ==
-  {
-    parallel::distributed::Vector<double> v2 (v);
-    bool equal = (v2 == v);
-    if (myid == 0)
-      deallog << " v==v2 ? " << equal << std::endl;
-
-    bool not_equal = (v2 != v);
-    if (myid == 0)
-      deallog << " v!=v2 ? " << not_equal << std::endl;
-
-    // change v2 on one proc only
-    if (myid == 0)
-      v2.local_element(1) = 2.2212;
-
-    equal = (v2 == v);
-    if (myid == 0)
-      deallog << " v==v2 ? " << equal << std::endl;
-    not_equal = (v2 != v);
-    if (myid == 0)
-      deallog << " v!=v2 ? " << not_equal << std::endl;
-
-    // reset
-    v2 = v;
-    equal = (v2 == v);
-    if (myid == 0)
-      deallog << " v==v2 ? " << equal << std::endl;
-    not_equal = (v2 != v);
-    if (myid == 0)
-      deallog << " v!=v2 ? " << not_equal << std::endl;
-
-    // change some value on all procs
-    if (myid < 8)
-      v2.local_element(0) = -1;
-    equal = (v2 == v);
-    if (myid == 0)
-      deallog << " v==v2 ? " << equal << std::endl;
-    not_equal = (v2 != v);
-    if (myid == 0)
-      deallog << " v!=v2 ? " << not_equal << std::endl;
-  }
-
   // check all_zero
   {
     bool allzero = v.all_zero();
@@ -181,38 +139,6 @@ void test ()
       deallog << " v2==0 ? " << allzero << std::endl;
   }
 
-
-  // check all_non_negative
-  {
-    bool allnonneg = v.is_non_negative();
-    if (myid == 0)
-      deallog << " v>=0 ? " << allnonneg << std::endl;
-    parallel::distributed::Vector<double> v2, v3;
-
-    // vector where all processors have
-    // non-negative entries
-    v2 = v;
-    if (myid < 8)
-      v2.local_element(0) = -1;
-    allnonneg = v2.is_non_negative();
-    if (myid == 0)
-      deallog << " v2>=0 ? " << allnonneg << std::endl;
-
-    // zero vector
-    v3.reinit (v2);
-    allnonneg = v3.is_non_negative();
-    if (myid == 0)
-      deallog << " v3>=0 ? " << allnonneg << std::endl;
-
-    // only one processor has non-negative entry
-    v3 = v;
-    if (myid == 1 || numproc==1)
-      v3.local_element(0) = -1;
-    allnonneg = v3.is_non_negative();
-    if (myid == 0)
-      deallog << " v3>=0 ? " << allnonneg << std::endl;
-  }
-
   if (myid == 0)
     deallog << "OK" << std::endl;
 }
index 8ea2068ba255d16afa6f12316e4dc3708eec42e6..b86e4d10e685243636eb3defd0e8379ad797b84b 100644 (file)
@@ -6,19 +6,7 @@ DEAL:0::linfty norm: 2.000
 DEAL:0::l2.2 norm: 2.000
 DEAL:0::Mean value: 1.000
 DEAL:0::Inner product: 4.000
-DEAL:0:: v==v2 ? 1
-DEAL:0:: v!=v2 ? 0
-DEAL:0:: v==v2 ? 0
-DEAL:0:: v!=v2 ? 1
-DEAL:0:: v==v2 ? 1
-DEAL:0:: v!=v2 ? 0
-DEAL:0:: v==v2 ? 0
-DEAL:0:: v!=v2 ? 1
 DEAL:0:: v==0 ? 0
 DEAL:0:: v2==0 ? 1
 DEAL:0:: v2==0 ? 0
-DEAL:0:: v>=0 ? 1
-DEAL:0:: v2>=0 ? 0
-DEAL:0:: v3>=0 ? 1
-DEAL:0:: v3>=0 ? 0
 DEAL:0::OK
index d147bb71a6b7291bef7119401cdd32b6f95cd1db..9684a3ed6ea75e50d3fd6c8da16851c1917dffc0 100644 (file)
@@ -6,19 +6,7 @@ DEAL:0::linfty norm: 30.00
 DEAL:0::l2.2 norm: 63.47
 DEAL:0::Mean value: 15.00
 DEAL:0::Inner product: 2608.
-DEAL:0:: v==v2 ? 1
-DEAL:0:: v!=v2 ? 0
-DEAL:0:: v==v2 ? 0
-DEAL:0:: v!=v2 ? 1
-DEAL:0:: v==v2 ? 1
-DEAL:0:: v!=v2 ? 0
-DEAL:0:: v==v2 ? 0
-DEAL:0:: v!=v2 ? 1
 DEAL:0:: v==0 ? 0
 DEAL:0:: v2==0 ? 1
 DEAL:0:: v2==0 ? 0
-DEAL:0:: v>=0 ? 1
-DEAL:0:: v2>=0 ? 0
-DEAL:0:: v3>=0 ? 1
-DEAL:0:: v3>=0 ? 0
 DEAL:0::OK
index d16831617eb7f310f82d902588eb49967bd40050..2c6ee43693c5cfad9241dac11cd111b3577a891d 100644 (file)
@@ -6,19 +6,7 @@ DEAL:0::linfty norm: 14.00
 DEAL:0::l2.2 norm: 22.04
 DEAL:0::Mean value: 7.000
 DEAL:0::Inner product: 312.0
-DEAL:0:: v==v2 ? 1
-DEAL:0:: v!=v2 ? 0
-DEAL:0:: v==v2 ? 0
-DEAL:0:: v!=v2 ? 1
-DEAL:0:: v==v2 ? 1
-DEAL:0:: v!=v2 ? 0
-DEAL:0:: v==v2 ? 0
-DEAL:0:: v!=v2 ? 1
 DEAL:0:: v==0 ? 0
 DEAL:0:: v2==0 ? 1
 DEAL:0:: v2==0 ? 0
-DEAL:0:: v>=0 ? 1
-DEAL:0:: v2>=0 ? 0
-DEAL:0:: v3>=0 ? 1
-DEAL:0:: v3>=0 ? 0
 DEAL:0::OK
index 97865e55403219fc966236b958409bb705231733..931cc4926cc636baa64cb4a971659f8bf6985f1c 100644 (file)
@@ -68,7 +68,7 @@ void test ()
   for (unsigned i=0; i<local_size; ++i)
     w.local_element(i) = 2.0 * (i + my_start);
 
-  v.copy_from(w);
+  v = w;
   v.update_ghost_values();
 
   // check local values for correctness
@@ -96,10 +96,11 @@ void test ()
       AssertThrow (v(ghost_indices[i])==v.local_element(local_size+i-5), ExcInternalError());
 
 
-  // now the same again, but import ghosts
-  // through the call to copy_from
+  // now the same again, but import ghosts automatically because v had ghosts
+  // set before calling operator =
   v.reinit (local_owned, local_relevant, MPI_COMM_WORLD);
-  v.copy_from(w, true);
+  v.update_ghost_values();
+  v = w;
 
   // check local values for correctness
   for (unsigned int i=0; i<local_size; ++i)
@@ -125,17 +126,6 @@ void test ()
     for (unsigned int i=5; i<10; ++i)
       AssertThrow (v(ghost_indices[i])==v.local_element(local_size+i-5), ExcInternalError());
 
-
-  // now do not call import_ghosts and check
-  // whether ghosts really are zero
-  v.reinit (local_owned, local_relevant, MPI_COMM_WORLD);
-  v.copy_from(w, false);
-
-  // check non-local entries on all processors
-  for (unsigned int i=0; i<10; ++i)
-    if (local_owned.is_element (ghost_indices[i]) == false)
-      AssertThrow (v(ghost_indices[i]) == 0., ExcInternalError());
-
   if (myid == 0)
     deallog << "OK" << std::endl;
 }
index ead41867912752fdcd6da19a674c84d572e92595..e9bc30d53e8276572cd7bab3bf46eded3d39ec6a 100644 (file)
@@ -73,7 +73,7 @@ void test ()
   if (myid==0) deallog << "OK" << std::endl;
 
   if (myid==0) deallog << "Check add (vector): ";
-  y.add (w);
+  y.add (1., w);
   for (int i=0; i<actual_local_size; ++i)
     AssertThrow (y.local_element(i) == 3*(i+my_start)+1042, ExcInternalError());
   if (myid==0) deallog << "OK" << std::endl;
@@ -104,7 +104,8 @@ void test ()
   if (myid==0) deallog << "OK" << std::endl;
 
   if (myid==0) deallog << "Check sadd (factor, factor, vector, factor, vector, factor, vector): ";
-  y.sadd (-1.,1.,v, 2., w, 2., x);
+  y.sadd (-1.,1.,v, 2., w);
+  y.add(2., x);
   for (int i=0; i<actual_local_size; ++i)
     AssertThrow (y.local_element(i) == 20000, ExcInternalError());
   if (myid==0) deallog << "OK" << std::endl;
@@ -147,7 +148,8 @@ void test ()
   if (myid==0) deallog << "OK" << std::endl;
 
   if (myid==0) deallog << "Check equ (factor, vector, factor, vector, factor, vector): ";
-  y. equ (10., v, -2., w, 3., x);
+  y. equ (10., v, -2., w);
+  y. add (3., x);
   for (int i=0; i<actual_local_size; ++i)
     AssertThrow (y.local_element(i) == 6.*(i+my_start)+28000, ExcInternalError());
   if (myid==0) deallog << "OK" << std::endl;
@@ -155,7 +157,7 @@ void test ()
   if (myid==0) deallog << "Check equ<float> (factor, vector): ";
   parallel::distributed::Vector<float> z;
   z = v;
-  y.equ (1., z);
+  y = z;
   for (int i=0; i<actual_local_size; ++i)
     AssertThrow (y.local_element(i) == i+my_start, ExcInternalError());
   if (myid==0) deallog << "OK" << std::endl;

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.