]> https://gitweb.dealii.org/ - dealii.git/commitdiff
bug fix applyBinary function in Rol::VectorAdaptor 5702/head
authorvishalkenchan <vishal.boddu@fau.de>
Fri, 5 Jan 2018 13:40:06 +0000 (14:40 +0100)
committervishalkenchan <vishal.boddu@fau.de>
Fri, 5 Jan 2018 13:40:06 +0000 (14:40 +0100)
include/deal.II/optimization/rol/vector_adaptor.h
tests/rol/vector_adaptor_apply_binary_01.cc [new file with mode: 0644]
tests/rol/vector_adaptor_apply_binary_01.output [new file with mode: 0644]

index cc068a2ba5f42ff9e22bbb57905bfa11875c1f9a..6264742d79f6d4b33e1ce6b5007d2c4abf257cba 100644 (file)
@@ -229,8 +229,8 @@ namespace Rol
      * Apply binary function @p f along with ROL::Vector @p rol_vector to all
      * the elements of the wrapped vector.
      */
-    void applyBinary (const ROL::Elementwise::UnaryFunction<value_type> &f,
-                      const ROL::Vector<value_type>                     &rol_vector);
+    void applyBinary (const ROL::Elementwise::BinaryFunction<value_type> &f,
+                      const ROL::Vector<value_type>                      &rol_vector);
 
     /**
      * Return the accumulated value on applying reduction operation @p r on
@@ -420,8 +420,8 @@ namespace Rol
   template<typename VectorType>
   void
   VectorAdaptor<VectorType>::
-  applyBinary (const ROL::Elementwise::UnaryFunction<value_type> &f,
-               const ROL::Vector<value_type>                     &rol_vector)
+  applyBinary (const ROL::Elementwise::BinaryFunction<value_type> &f,
+               const ROL::Vector<value_type>                      &rol_vector)
   {
     Assert (this->dimension() == rol_vector.dimension(),
             ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
@@ -431,14 +431,13 @@ namespace Rol
 
     const VectorType &given_rol_vector = *(vector_adaptor.getVector());
 
-    const typename VectorType::iterator
-    vend   = vector_ptr->end(),
-    rolend = given_rol_vector.end();
+    const typename VectorType::iterator       vend   = vector_ptr->end();
+    const typename VectorType::const_iterator rolend = given_rol_vector.end();
 
-    for (typename VectorType::iterator
-         l_iterator  = vector_ptr->begin(), r_iterator  = given_rol_vector.begin();
-         l_iterator != vend              && r_iterator != rolend;
-         l_iterator++,                      r_iterator++)
+    typename VectorType::const_iterator r_iterator = given_rol_vector.begin();
+    for (typename VectorType::iterator  l_iterator = vector_ptr->begin();
+         l_iterator != vend && r_iterator != rolend;
+         l_iterator++,         r_iterator++)
       *l_iterator = f.apply(*l_iterator, *r_iterator);
 
     vector_ptr->compress (VectorOperation::insert);
diff --git a/tests/rol/vector_adaptor_apply_binary_01.cc b/tests/rol/vector_adaptor_apply_binary_01.cc
new file mode 100644 (file)
index 0000000..769b26d
--- /dev/null
@@ -0,0 +1,150 @@
+// ---------------------------------------------------------------------
+//
+// 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 the Rol::VectorAdaptor's applyBinary function applied to a fully
+// distributed (MPI) vector.
+
+#include "../tests.h"
+
+#include <deal.II/lac/generic_linear_algebra.h>
+#include <deal.II/optimization/rol/vector_adaptor.h>
+
+using namespace dealii;
+
+// Taken from deal.II's test: parallel_vector_07
+template <typename VectorType>
+void prepare_vector (VectorType &v)
+{
+  const unsigned int
+  myid    = dealii::Utilities::MPI::this_mpi_process (MPI_COMM_WORLD),
+  numproc = dealii::Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD);
+
+  const unsigned int set = 10;
+  AssertIndexRange (numproc, set-2);
+  const unsigned int local_size = set - myid;
+  unsigned int global_size = 0;
+  unsigned int my_start = 0;
+  for (unsigned int i=0; i<numproc; ++i)
+    {
+      global_size += set - i;
+      if (i<myid)
+        my_start += set - i;
+    }
+  // each processor owns some indices and all
+  // are ghosting elements from three
+  // processors (the second). some entries
+  // are right around the border between two
+  // processors
+  IndexSet local_owned(global_size);
+  local_owned.add_range(my_start, my_start + local_size);
+
+  // --- Prepare vector.
+  v.reinit (local_owned, MPI_COMM_WORLD);
+}
+
+
+template <typename VectorType>
+void test ()
+{
+  VectorType a, b;
+  prepare_vector (a);
+  prepare_vector (b);
+
+  Testing::srand(1);
+
+  for (auto iterator = a.begin(); iterator != a.end(); iterator++)
+    *iterator = static_cast<double>(Testing::rand())
+                /
+                RAND_MAX;
+
+  for (auto iterator = b.begin(); iterator != b.end(); iterator++)
+    *iterator = static_cast<double>(Testing::rand())
+                /
+                RAND_MAX;
+
+  a.compress(VectorOperation::insert);
+  b.compress(VectorOperation::insert);
+
+  Teuchos::RCP<VectorType> a_rcp (new VectorType(a));
+  Teuchos::RCP<VectorType> b_rcp (new VectorType(b));
+
+  ROL::Elementwise::Multiply<double> mult;
+  ROL::Elementwise::Plus<double>     plus;
+
+  // --- Testing the constructor
+  Rol::VectorAdaptor<VectorType> a_rol (a_rcp);
+  Rol::VectorAdaptor<VectorType> b_rol (b_rcp);
+
+  a_rol.print(std::cout);
+  b_rol.print(std::cout);
+
+  a_rol.applyBinary(mult, b_rol);
+  a_rol.print(std::cout);
+
+  a_rol.applyBinary(plus, b_rol);
+  a_rol.print(std::cout);
+}
+
+
+
+int main (int argc, char **argv)
+{
+
+
+  dealii::Utilities::MPI::MPI_InitFinalize
+  mpi_initialization (argc,
+                      argv,
+                      1);
+
+  unsigned int myid = dealii::Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+  deallog.push(dealii::Utilities::int_to_string(myid));
+
+
+  if (myid == 0)
+    {
+      deallog.depth_console(10); // initlog();
+      deallog << std::setprecision(4);
+    }
+
+  try
+    {
+      test<LinearAlgebraTrilinos::MPI::Vector>();
+      test<LinearAlgebra::distributed::Vector<double>>();
+    }
+  catch (std::exception &exc)
+    {
+      std::cerr << std::endl << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+      std::cerr << "Exception on processing: " << std::endl
+                << exc.what() << std::endl
+                << "Aborting!" << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+
+      return 1;
+    }
+  catch (...)
+    {
+      std::cerr << std::endl << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+      std::cerr << "Unknown exception!" << std::endl
+                << "Aborting!" << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+      return 1;
+    };
+}
diff --git a/tests/rol/vector_adaptor_apply_binary_01.output b/tests/rol/vector_adaptor_apply_binary_01.output
new file mode 100644 (file)
index 0000000..c903498
--- /dev/null
@@ -0,0 +1,20 @@
+8.402e-01 3.944e-01 7.831e-01 7.984e-01 9.116e-01 1.976e-01 3.352e-01 7.682e-01 2.778e-01 5.540e-01 
+4.774e-01 6.289e-01 3.648e-01 5.134e-01 9.522e-01 9.162e-01 6.357e-01 7.173e-01 1.416e-01 6.070e-01 
+4.011e-01 2.480e-01 2.857e-01 4.099e-01 8.681e-01 1.810e-01 2.131e-01 5.510e-01 3.933e-02 3.362e-01 
+8.785e-01 8.769e-01 6.504e-01 9.233e-01 1.820e+00 1.097e+00 8.488e-01 1.268e+00 1.809e-01 9.432e-01 
+Process #0
+Local range: [0, 10), global size: 10
+Vector data:
+8.402e-01 3.944e-01 7.831e-01 7.984e-01 9.116e-01 1.976e-01 3.352e-01 7.682e-01 2.778e-01 5.540e-01 
+Process #0
+Local range: [0, 10), global size: 10
+Vector data:
+4.774e-01 6.289e-01 3.648e-01 5.134e-01 9.522e-01 9.162e-01 6.357e-01 7.173e-01 1.416e-01 6.070e-01 
+Process #0
+Local range: [0, 10), global size: 10
+Vector data:
+4.011e-01 2.480e-01 2.857e-01 4.099e-01 8.681e-01 1.810e-01 2.131e-01 5.510e-01 3.933e-02 3.362e-01 
+Process #0
+Local range: [0, 10), global size: 10
+Vector data:
+8.785e-01 8.769e-01 6.504e-01 9.233e-01 1.820e+00 1.097e+00 8.488e-01 1.268e+00 1.809e-01 9.432e-01 

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.