]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add an 'add' method to the native parallel vectors.
authorDavid Wells <wellsd2@rpi.edu>
Wed, 28 Dec 2016 17:39:22 +0000 (12:39 -0500)
committerDavid Wells <wellsd2@rpi.edu>
Mon, 9 Jan 2017 14:40:24 +0000 (09:40 -0500)
This is necessary for an upcoming patch that uses the method inside
ConstraintMatrix::distribute_local_to_global.

12 files changed:
include/deal.II/lac/la_parallel_block_vector.h
include/deal.II/lac/la_parallel_block_vector.templates.h
include/deal.II/lac/la_parallel_vector.h
include/deal.II/lac/la_parallel_vector.templates.h
tests/mpi/parallel_block_vector_04.cc [new file with mode: 0644]
tests/mpi/parallel_block_vector_04.mpirun=1.output [new file with mode: 0644]
tests/mpi/parallel_block_vector_04.mpirun=10.output [new file with mode: 0644]
tests/mpi/parallel_block_vector_04.mpirun=4.output [new file with mode: 0644]
tests/mpi/parallel_vector_19.cc [new file with mode: 0644]
tests/mpi/parallel_vector_19.mpirun=1.output [new file with mode: 0644]
tests/mpi/parallel_vector_19.mpirun=10.output [new file with mode: 0644]
tests/mpi/parallel_vector_19.mpirun=4.output [new file with mode: 0644]

index f96e64052e503700e326fe6517f8ea9c4adf74ef..bbfeae8134e2337c579dd320c18f2b72f2719a52 100644 (file)
@@ -458,6 +458,13 @@ namespace LinearAlgebra
       virtual void add(const Number a, const VectorSpaceVector<Number> &V,
                        const Number b, const VectorSpaceVector<Number> &W);
 
+      /**
+       * A collective add operation: This function adds a whole set of values
+       * stored in @p values to the vector components specified by @p indices.
+       */
+      virtual void add (const std::vector<size_type> &indices,
+                        const std::vector<Number>    &values);
+
       /**
        * Scaling and simple addition of a multiple of a vector, i.e. <tt>*this =
        * s*(*this)+a*V</tt>.
index 1f5ea732de6c7077dce7ba7ac8f7a38da1564e4f..116c34d2050af0c2285883c016be9443a5ba1332 100644 (file)
@@ -497,6 +497,18 @@ namespace LinearAlgebra
 
 
 
+    template <typename Number>
+    void
+    BlockVector<Number>::add (const std::vector<size_type> &indices,
+                              const std::vector<Number>    &values)
+    {
+      for (size_type i=0; i<indices.size(); ++i)
+        (*this)(indices[i]) += values[i];
+    }
+
+
+
+
     template <typename Number>
     bool
     BlockVector<Number>::all_zero () const
index 43fa05a2641011adbe36c27fe1cefbd9258456b7..b73fd4ea7ec61e588bebe01f381785ec96a14c67 100644 (file)
@@ -584,6 +584,13 @@ namespace LinearAlgebra
       virtual void add(const Number a, const VectorSpaceVector<Number> &V,
                        const Number b, const VectorSpaceVector<Number> &W);
 
+      /**
+       * A collective add operation: This function adds a whole set of values
+       * stored in @p values to the vector components specified by @p indices.
+       */
+      virtual void add (const std::vector<size_type> &indices,
+                        const std::vector<Number>    &values);
+
       /**
        * Scaling and simple addition of a multiple of a vector, i.e. <tt>*this =
        * s*(*this)+a*V</tt>.
index 37d8433162b04cdcd7cbd220beeaff493a55dd62..48a24d59da19cc71a0443347bf5c70b814eaa218 100644 (file)
@@ -1044,6 +1044,20 @@ namespace LinearAlgebra
 
 
 
+    template <typename Number>
+    void
+    Vector<Number>::add (const std::vector<size_type> &indices,
+                         const std::vector<Number>    &values)
+    {
+      for (std::size_t i=0; i<indices.size(); ++i)
+        {
+          this->operator()(indices[i]) += values[i];
+        }
+    }
+
+
+
+
     template <typename Number>
     void
     Vector<Number>::sadd (const Number x,
diff --git a/tests/mpi/parallel_block_vector_04.cc b/tests/mpi/parallel_block_vector_04.cc
new file mode 100644 (file)
index 0000000..d6eb2f5
--- /dev/null
@@ -0,0 +1,200 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 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.
+//
+// ---------------------------------------------------------------------
+
+
+// Check global reduction operations (norms, operator==, operator!=) with the
+// parallel block vector class. This is the same as parallel_vector_01 except
+// that this test uses the version of the add() method that takes a std::vector
+// of indices and a std::vector of values to build the block vector instead of
+// vector assignment.
+
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+#include <deal.II/base/index_set.h>
+#include <deal.II/lac/la_parallel_block_vector.h>
+#include <deal.II/lac/la_parallel_vector.h>
+#include <fstream>
+#include <iostream>
+#include <vector>
+
+
+void test ()
+{
+  unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+  unsigned int numproc = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD);
+
+  if (myid==0) deallog << "numproc=" << numproc << std::endl;
+
+
+  // each processor from processor 1 to 8
+  // owns 2 indices (the other processors do
+  // not own any dof), and all processors are
+  // ghosting element 1 (the second)
+  IndexSet local_owned(std::min(16U,numproc*2));
+  if (myid < 8)
+    local_owned.add_range(myid*2,myid*2+2);
+  IndexSet local_relevant(numproc*2);
+  local_relevant = local_owned;
+  local_relevant.add_range(1,2);
+
+  LinearAlgebra::distributed::Vector<double> v(local_owned, local_relevant,
+                                               MPI_COMM_WORLD);
+  std::vector<types::global_dof_index> indices;
+  std::vector<double> values;
+
+  // set local values
+  if (myid < 8)
+    {
+      // unlike the first test, we scale values here since these are added to
+      // w directly
+      indices.push_back(myid*2);
+      values.push_back(myid*4.0);
+      indices.push_back(myid*2 + 1);
+      values.push_back(myid*4.0 + 2.0);
+    }
+
+  LinearAlgebra::distributed::BlockVector<double> w(3);
+  for (unsigned int i=0; i<3; ++i)
+    w.block(i) = v;
+  w.collect_sizes();
+
+  // set up v for comparison purposes below
+  v.add(indices, values);
+  v.compress(VectorOperation::add);
+  // in analogy to the first test, v is scaled differently than w
+
+  if (myid < 8)
+    {
+      for (unsigned int i = 0; i<3; ++i)
+        {
+          // update the values in w
+          w.add(indices, values);
+          // update indices to point to the next block
+          indices[0] += v.size();
+          indices[1] += v.size();
+        }
+    }
+  w.compress(VectorOperation::add);
+
+  // check l2 norm
+  {
+    const double l2_norm = w.l2_norm();
+    if (myid == 0)
+      deallog << "l2 norm: " << l2_norm << std::endl;
+    AssertThrow (std::abs(v.l2_norm()*std::sqrt(3.)-w.l2_norm()) < 1e-13,
+                 ExcInternalError());
+  }
+
+  // check l1 norm
+  {
+    const double l1_norm = w.l1_norm();
+    if (myid == 0)
+      deallog << "l1 norm: " << l1_norm << std::endl;
+    AssertThrow (std::abs(v.l1_norm()*3.-w.l1_norm()) < 1e-14,
+                 ExcInternalError());
+  }
+
+  // check linfty norm
+  {
+    const double linfty_norm = w.linfty_norm();
+    if (myid == 0)
+      deallog << "linfty norm: " << linfty_norm << std::endl;
+    AssertThrow (v.linfty_norm()==w.linfty_norm(),
+                 ExcInternalError());
+  }
+
+  // check lp norm
+  {
+    const double lp_norm = w.lp_norm(2.2);
+    if (myid == 0)
+      deallog << "l2.2 norm: " << lp_norm << std::endl;
+
+    AssertThrow (std::fabs (w.l2_norm() - w.lp_norm(2.0)) < 1e-14,
+                 ExcInternalError());
+  }
+
+  // check mean value (should be equal to l1
+  // norm divided by vector size here since we
+  // have no negative entries)
+  {
+    const double mean = w.mean_value();
+    if (myid == 0)
+      deallog << "Mean value: " << mean << std::endl;
+
+    AssertThrow (std::fabs (mean * w.size() - w.l1_norm()) < 1e-15,
+                 ExcInternalError());
+  }
+  // check inner product
+  {
+    const double norm_sqr = w.l2_norm() * w.l2_norm();
+    AssertThrow (std::fabs(w * w - norm_sqr) < 1e-12,
+                 ExcInternalError());
+    LinearAlgebra::distributed::BlockVector<double> w2;
+    w2 = w;
+    AssertThrow (std::fabs(w2 * w - norm_sqr) < 1e-12,
+                 ExcInternalError());
+
+    if (myid<8)
+      w2.block(0).local_element(0) = -1;
+    const double inner_prod = w * w2;
+    if (myid == 0)
+      deallog << "Inner product: " << inner_prod << std::endl;
+  }
+
+  // check all_zero
+  {
+    bool allzero = w.all_zero();
+    if (myid == 0)
+      deallog << " v==0 ? " << allzero << std::endl;
+    LinearAlgebra::distributed::BlockVector<double> w2;
+    w2.reinit (w);
+    allzero = w2.all_zero();
+    if (myid == 0)
+      deallog << " v2==0 ? " << allzero << std::endl;
+
+    // now change one element to nonzero
+    if (myid == 0)
+      w2.block(1).local_element(1) = 1;
+    allzero = w2.all_zero();
+    if (myid == 0)
+      deallog << " v2==0 ? " << allzero << std::endl;
+  }
+
+  if (myid == 0)
+    deallog << "OK" << std::endl;
+}
+
+
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads());
+
+  unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+  deallog.push(Utilities::int_to_string(myid));
+
+  if (myid == 0)
+    {
+      std::ofstream logfile("output");
+      deallog.attach(logfile);
+      deallog << std::setprecision(4);
+      deallog.threshold_double(1.e-10);
+
+      test();
+    }
+  else
+    test();
+
+}
diff --git a/tests/mpi/parallel_block_vector_04.mpirun=1.output b/tests/mpi/parallel_block_vector_04.mpirun=1.output
new file mode 100644 (file)
index 0000000..3a3e518
--- /dev/null
@@ -0,0 +1,12 @@
+
+DEAL:0::numproc=1
+DEAL:0::l2 norm: 3.464
+DEAL:0::l1 norm: 6.000
+DEAL:0::linfty norm: 2.000
+DEAL:0::l2.2 norm: 3.295
+DEAL:0::Mean value: 1.000
+DEAL:0::Inner product: 12.00
+DEAL:0:: v==0 ? 0
+DEAL:0:: v2==0 ? 1
+DEAL:0:: v2==0 ? 0
+DEAL:0::OK
diff --git a/tests/mpi/parallel_block_vector_04.mpirun=10.output b/tests/mpi/parallel_block_vector_04.mpirun=10.output
new file mode 100644 (file)
index 0000000..ba129a4
--- /dev/null
@@ -0,0 +1,12 @@
+
+DEAL:0::numproc=10
+DEAL:0::l2 norm: 122.0
+DEAL:0::l1 norm: 720.0
+DEAL:0::linfty norm: 30.00
+DEAL:0::l2.2 norm: 104.6
+DEAL:0::Mean value: 15.00
+DEAL:0::Inner product: 1.253e+04
+DEAL:0:: v==0 ? 0
+DEAL:0:: v2==0 ? 1
+DEAL:0:: v2==0 ? 0
+DEAL:0::OK
diff --git a/tests/mpi/parallel_block_vector_04.mpirun=4.output b/tests/mpi/parallel_block_vector_04.mpirun=4.output
new file mode 100644 (file)
index 0000000..9ffa136
--- /dev/null
@@ -0,0 +1,12 @@
+
+DEAL:0::numproc=4
+DEAL:0::l2 norm: 40.99
+DEAL:0::l1 norm: 168.0
+DEAL:0::linfty norm: 14.00
+DEAL:0::l2.2 norm: 36.31
+DEAL:0::Mean value: 7.000
+DEAL:0::Inner product: 1432.
+DEAL:0:: v==0 ? 0
+DEAL:0:: v2==0 ? 1
+DEAL:0:: v2==0 ? 0
+DEAL:0::OK
diff --git a/tests/mpi/parallel_vector_19.cc b/tests/mpi/parallel_vector_19.cc
new file mode 100644 (file)
index 0000000..fcd5f0e
--- /dev/null
@@ -0,0 +1,96 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+
+// use same partition as in parallel_vector_06 to check the version of the
+// add() method that takes a std::vector of indices and a std::vector of
+// values.
+
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+#include <deal.II/base/index_set.h>
+#include <deal.II/lac/la_parallel_vector.h>
+#include <fstream>
+#include <iostream>
+#include <vector>
+
+
+void test ()
+{
+  unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+  unsigned int numproc = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD);
+
+  if (myid==0) deallog << "numproc=" << numproc << std::endl;
+
+
+  // each processor from processor 1 to 8
+  // owns 2 indices (the other processors do
+  // not own any dof), and all processors are
+  // ghosting element 1 (the second)
+  IndexSet local_owned(std::min(16U,numproc*2));
+  if (myid < 8)
+    local_owned.add_range(myid*2,myid*2+2);
+  IndexSet local_relevant(numproc*2);
+  local_relevant = local_owned;
+  local_relevant.add_range(1,2);
+
+  LinearAlgebra::distributed::Vector<double> v(local_owned, local_owned, MPI_COMM_WORLD);
+
+  // set local values
+  if (myid < 8)
+    {
+      types::global_dof_index n_elements = 2;
+      std::vector<types::global_dof_index> indices(n_elements);
+      indices[0] = myid*2;
+      indices[1] = myid*2+1;
+      std::vector<double> values(n_elements);
+      values[0] = myid*2.0;
+      values[1] = myid*2.0+1.0;
+      v.add (indices, values);
+    }
+  v.compress(VectorOperation::insert);
+  v*=2.0;
+  if (myid < 8)
+    {
+      AssertThrow(v(myid*2) == myid*4.0, ExcInternalError());
+      AssertThrow(v(myid*2+1) == myid*4.0+2.0, ExcInternalError());
+    }
+
+  if (myid == 0)
+    deallog << "OK" << std::endl;
+}
+
+
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads());
+
+  unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+  deallog.push(Utilities::int_to_string(myid));
+
+  if (myid == 0)
+    {
+      std::ofstream logfile("output");
+      deallog.attach(logfile);
+      deallog << std::setprecision(4);
+      deallog.threshold_double(1.e-10);
+
+      test();
+    }
+  else
+    test();
+
+}
diff --git a/tests/mpi/parallel_vector_19.mpirun=1.output b/tests/mpi/parallel_vector_19.mpirun=1.output
new file mode 100644 (file)
index 0000000..bb1593d
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL:0::numproc=1
+DEAL:0::OK
diff --git a/tests/mpi/parallel_vector_19.mpirun=10.output b/tests/mpi/parallel_vector_19.mpirun=10.output
new file mode 100644 (file)
index 0000000..999baf1
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL:0::numproc=10
+DEAL:0::OK
diff --git a/tests/mpi/parallel_vector_19.mpirun=4.output b/tests/mpi/parallel_vector_19.mpirun=4.output
new file mode 100644 (file)
index 0000000..6f5d277
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL:0::numproc=4
+DEAL:0::OK

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.