]> https://gitweb.dealii.org/ - dealii.git/commitdiff
added rol tests 5380/head
authorvishalkenchan <vishal.boddu@fau.de>
Fri, 10 Nov 2017 09:15:36 +0000 (10:15 +0100)
committervishalkenchan <vishal.boddu@fau.de>
Fri, 10 Nov 2017 12:28:07 +0000 (13:28 +0100)
tests/rol/CMakeLists.txt [new file with mode: 0644]
tests/rol/vector_adaptor_01.cc [new file with mode: 0644]
tests/rol/vector_adaptor_01.output [new file with mode: 0644]
tests/rol/vector_adaptor_02.cc [new file with mode: 0644]
tests/rol/vector_adaptor_02.output [new file with mode: 0644]
tests/rol/vector_adaptor_no_ghost_01.cc [new file with mode: 0644]
tests/rol/vector_adaptor_no_ghost_01.mpirun=5.output [new file with mode: 0644]
tests/rol/vector_adaptor_no_ghost_01.output [new file with mode: 0644]
tests/rol/vector_adaptor_with_ghost_01.cc [new file with mode: 0644]
tests/rol/vector_adaptor_with_ghost_01.mpirun=10.output [new file with mode: 0644]
tests/rol/vector_adaptor_with_ghost_01.mpirun=7.output [new file with mode: 0644]

diff --git a/tests/rol/CMakeLists.txt b/tests/rol/CMakeLists.txt
new file mode 100644 (file)
index 0000000..57f73cd
--- /dev/null
@@ -0,0 +1,6 @@
+CMAKE_MINIMUM_REQUIRED(VERSION 2.8.12)
+INCLUDE(../setup_testsubproject.cmake)
+PROJECT(testsuite CXX)
+IF(DEAL_II_WITH_TRILINOS AND DEAL_II_TRILINOS_WITH_ROL)
+  DEAL_II_PICKUP_TESTS()
+ENDIF()
diff --git a/tests/rol/vector_adaptor_01.cc b/tests/rol/vector_adaptor_01.cc
new file mode 100644 (file)
index 0000000..45be97c
--- /dev/null
@@ -0,0 +1,105 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+// Test to check Rol::VectorAdaptor::set() and Rol::VectorAdaptor::plus().
+
+#include "../tests.h"
+
+#include <deal.II/lac/generic_linear_algebra.h>
+#include <deal.II/optimization/rol/vector_adaptor.h>
+
+using namespace dealii;
+
+template <typename VectorType>
+void test (const VectorType &given_vector)
+{
+  Teuchos::RCP<VectorType> given_vector_rcp (new VectorType(given_vector));
+
+  // --- Testing the constructor
+  Rol::VectorAdaptor<VectorType> given_vector_rol (given_vector_rcp);
+  AssertThrow (given_vector == *given_vector_rol.getVector(), ExcInternalError());
+
+
+  Teuchos::RCP<VectorType> w_rcp =  Teuchos::rcp (new VectorType);
+  Rol::VectorAdaptor<VectorType> w_rol (w_rcp);
+
+  // --- Testing VectorAdaptor::set()
+  {
+    w_rol.set(given_vector_rol);
+    AssertThrow (given_vector == *w_rol.getVector(), ExcInternalError());
+  }
+
+  // --- Testing VectorAdaptor::plus()
+  {
+    VectorType u;
+    u = given_vector;
+    u *= 2.;
+    w_rol.plus (given_vector_rol);
+    AssertThrow (u == *w_rol.getVector(), ExcInternalError());
+  }
+
+  deallog << "OK" << std::endl;
+}
+
+
+
+int main (int argc, char **argv)
+{
+  deallog.depth_console(10);
+
+  dealii::Utilities::MPI::MPI_InitFinalize
+  mpi_initialization (argc,
+                      argv,
+                      dealii::numbers::invalid_unsigned_int);
+
+  try
+    {
+      {
+        LinearAlgebraTrilinos::MPI::Vector trilinos_vector;
+        trilinos_vector.reinit(complete_index_set(100), MPI_COMM_WORLD);
+
+        // set the first vector
+        for (unsigned int i=0; i<trilinos_vector.size(); ++i)
+          trilinos_vector(i) = i;
+
+        test (trilinos_vector);
+
+      }
+    }
+  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_01.output b/tests/rol/vector_adaptor_01.output
new file mode 100644 (file)
index 0000000..559e3f5
--- /dev/null
@@ -0,0 +1 @@
+DEAL::OK
diff --git a/tests/rol/vector_adaptor_02.cc b/tests/rol/vector_adaptor_02.cc
new file mode 100644 (file)
index 0000000..5e80e81
--- /dev/null
@@ -0,0 +1,141 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+#include <cmath>
+#include <iostream>
+#include <sstream>
+
+#include <deal.II/lac/generic_linear_algebra.h>
+#include <deal.II/optimization/rol/vector_adaptor.h>
+
+#include "ROL_Objective.hpp"
+#include "ROL_Algorithm.hpp"
+#include "ROL_LineSearchStep.hpp"
+#include "ROL_StatusTest.hpp"
+#include "Teuchos_GlobalMPISession.hpp"
+
+// Use ROL to minimize the objective function, f(x,y) = x^2 + y^2.
+
+using namespace dealii;
+
+using VectorType = typename dealii::Vector<double>;
+
+template<class Real=double, typename Xprim=Rol::VectorAdaptor<VectorType> >
+class QuadraticObjective : public ROL::Objective<Real>
+{
+
+private:
+
+  Teuchos::RCP<const VectorType>
+  get_rcp_to_VectorType (const ROL::Vector<Real> &x)
+  {
+    return (Teuchos::dyn_cast<const Xprim>(x)).getVector();
+  }
+
+  Teuchos::RCP<dealii::Vector<Real> >
+  get_rcp_to_VectorType (ROL::Vector<Real> &x)
+  {
+    return (Teuchos::dyn_cast<Xprim>(x)).getVector();
+  }
+
+public:
+
+  Real value (const ROL::Vector<Real> &x,
+              Real                    &tol)
+  {
+    Assert (x.dimension()==2,
+            ExcInternalError());
+
+    return x.dot(x);
+  }
+
+  void gradient (ROL::Vector<Real>       &g,
+                 const ROL::Vector<Real> &x,
+                 Real                    &tol)
+  {
+    Teuchos::RCP<const VectorType> xp = this->get_rcp_to_VectorType(x);
+    Teuchos::RCP<VectorType> gp = this->get_rcp_to_VectorType(g);
+
+    (*gp)[0] = 2. * (*xp)[0];
+    (*gp)[1] = 2. * (*xp)[1];
+  }
+
+};
+
+void test (const double x,
+           const double y)
+{
+  typedef double RealT;
+
+  QuadraticObjective<RealT> quad_objective;
+
+  Teuchos::RCP<std::ostream> outStream = Teuchos::rcp(&std::cout, false);
+  Teuchos::RCP<VectorType> x_rcp       = Teuchos::rcp (new VectorType);
+
+  x_rcp->reinit (2);
+
+  (*x_rcp)[0] = x;
+  (*x_rcp)[1] = y;
+
+  Rol::VectorAdaptor<VectorType> x_rol(x_rcp);
+
+  Teuchos::ParameterList parlist;
+  // Set parameters.
+  parlist.sublist("Secant").set("Use as Preconditioner", false);
+  // Define algorithm.
+  ROL::Algorithm<RealT> algo("Line Search", parlist);
+
+  // Run Algorithm
+  algo.run(x_rol, quad_objective, true, *outStream);
+
+  Teuchos::RCP<const VectorType> xg = x_rol.getVector();
+  std::cout << "The solution to minimization problem is: ";
+  std::cout << (*xg)[0] << " " << (*xg)[1] << std::endl;
+}
+
+int main (int argc, char **argv)
+{
+  try
+    {
+      test(  10,  -2);
+      test(-0.1, 0.1);
+      test( 9.1,-6.1);
+    }
+  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;
+      throw;
+    }
+  catch (...)
+    {
+      std::cerr << std::endl << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+      std::cerr << "Unknown exception!" << std::endl
+                << "Aborting!" << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+      throw;
+    }
+
+  return 0;
+}
diff --git a/tests/rol/vector_adaptor_02.output b/tests/rol/vector_adaptor_02.output
new file mode 100644 (file)
index 0000000..75350d0
--- /dev/null
@@ -0,0 +1,21 @@
+
+Quasi-Newton Method with Limited-Memory BFGS
+Line Search: Cubic Interpolation satisfying Strong Wolfe Conditions
+  iter  value          gnorm          snorm          #fval     #grad     ls_#fval  ls_#grad  
+  0     1.040000e+02   2.039608e+01   
+  1     0.000000e+00   0.000000e+00   1.019804e+01   4         2         2         0         
+The solution to minimization problem is: 0 0
+
+Quasi-Newton Method with Limited-Memory BFGS
+Line Search: Cubic Interpolation satisfying Strong Wolfe Conditions
+  iter  value          gnorm          snorm          #fval     #grad     ls_#fval  ls_#grad  
+  0     2.000000e-02   2.828427e-01   
+  1     0.000000e+00   0.000000e+00   1.414214e-01   4         2         2         0         
+The solution to minimization problem is: 0 0
+
+Quasi-Newton Method with Limited-Memory BFGS
+Line Search: Cubic Interpolation satisfying Strong Wolfe Conditions
+  iter  value          gnorm          snorm          #fval     #grad     ls_#fval  ls_#grad  
+  0     1.200200e+02   2.191073e+01   
+  1     0.000000e+00   0.000000e+00   1.095536e+01   4         2         2         0         
+The solution to minimization problem is: 0 0
diff --git a/tests/rol/vector_adaptor_no_ghost_01.cc b/tests/rol/vector_adaptor_no_ghost_01.cc
new file mode 100644 (file)
index 0000000..b52ad2f
--- /dev/null
@@ -0,0 +1,155 @@
+// ---------------------------------------------------------------------
+//
+// 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 with MPI fully distributed vectors
+// using ROL::Vector's checkVector method.
+
+#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 = 200;
+  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, c;
+  prepare_vector (a);
+  prepare_vector (b);
+  prepare_vector (c);
+
+  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;
+
+  for (auto iterator = c.begin(); iterator != c.end(); iterator++)
+    *iterator = static_cast<double>(Testing::rand())
+                /
+                RAND_MAX;
+
+  a.compress(VectorOperation::insert);
+  b.compress(VectorOperation::insert);
+  c.compress(VectorOperation::insert);
+
+  Teuchos::RCP<VectorType> a_rcp (new VectorType(a));
+  Teuchos::RCP<VectorType> b_rcp (new VectorType(b));
+  Teuchos::RCP<VectorType> c_rcp (new VectorType(c));
+
+  // --- Testing the constructor
+  Rol::VectorAdaptor<VectorType> a_rol (a_rcp);
+  Rol::VectorAdaptor<VectorType> b_rol (b_rcp);
+  Rol::VectorAdaptor<VectorType> c_rol (c_rcp);
+
+  Teuchos::RCP<std::ostream> out_stream;
+  Teuchos::oblackholestream bhs; // outputs nothing
+
+  if (dealii::Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0)
+    out_stream = Teuchos::rcp(&std::cout, false);
+  else
+    out_stream = Teuchos::rcp(&bhs, false);
+
+  a_rol.checkVector (b_rol, c_rol, true, *out_stream);
+}
+
+
+
+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_no_ghost_01.mpirun=5.output b/tests/rol/vector_adaptor_no_ghost_01.mpirun=5.output
new file mode 100644 (file)
index 0000000..6b11695
--- /dev/null
@@ -0,0 +1,36 @@
+
+********** Begin verification of linear algebra. *********************************************
+
+Commutativity of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Associativity of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Identity element of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Inverse elements of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Identity element of scalar multiplication. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Consistency of scalar multiplication with field multiplication. Consistency error: >>>>>>>>>>> 0
+Distributivity of scalar multiplication with respect to field addition. Consistency error: >>> 0
+Distributivity of scalar multiplication with respect to vector addition. Consistency error: >> 0
+Commutativity of dot (inner) product over the field of reals. Consistency error: >>>>>>>>>>>>> 0
+Additivity of dot (inner) product. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Consistency of scalar multiplication and norm. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Reflexivity. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+
+********** End verification of linear algebra. ***********************************************
+
+
+********** Begin verification of linear algebra. *********************************************
+
+Commutativity of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Associativity of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Identity element of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Inverse elements of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Identity element of scalar multiplication. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Consistency of scalar multiplication with field multiplication. Consistency error: >>>>>>>>>>> 0
+Distributivity of scalar multiplication with respect to field addition. Consistency error: >>> 0
+Distributivity of scalar multiplication with respect to vector addition. Consistency error: >> 0
+Commutativity of dot (inner) product over the field of reals. Consistency error: >>>>>>>>>>>>> 0
+Additivity of dot (inner) product. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Consistency of scalar multiplication and norm. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Reflexivity. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+
+********** End verification of linear algebra. ***********************************************
+
diff --git a/tests/rol/vector_adaptor_no_ghost_01.output b/tests/rol/vector_adaptor_no_ghost_01.output
new file mode 100644 (file)
index 0000000..6b11695
--- /dev/null
@@ -0,0 +1,36 @@
+
+********** Begin verification of linear algebra. *********************************************
+
+Commutativity of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Associativity of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Identity element of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Inverse elements of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Identity element of scalar multiplication. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Consistency of scalar multiplication with field multiplication. Consistency error: >>>>>>>>>>> 0
+Distributivity of scalar multiplication with respect to field addition. Consistency error: >>> 0
+Distributivity of scalar multiplication with respect to vector addition. Consistency error: >> 0
+Commutativity of dot (inner) product over the field of reals. Consistency error: >>>>>>>>>>>>> 0
+Additivity of dot (inner) product. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Consistency of scalar multiplication and norm. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Reflexivity. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+
+********** End verification of linear algebra. ***********************************************
+
+
+********** Begin verification of linear algebra. *********************************************
+
+Commutativity of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Associativity of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Identity element of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Inverse elements of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Identity element of scalar multiplication. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Consistency of scalar multiplication with field multiplication. Consistency error: >>>>>>>>>>> 0
+Distributivity of scalar multiplication with respect to field addition. Consistency error: >>> 0
+Distributivity of scalar multiplication with respect to vector addition. Consistency error: >> 0
+Commutativity of dot (inner) product over the field of reals. Consistency error: >>>>>>>>>>>>> 0
+Additivity of dot (inner) product. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Consistency of scalar multiplication and norm. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Reflexivity. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+
+********** End verification of linear algebra. ***********************************************
+
diff --git a/tests/rol/vector_adaptor_with_ghost_01.cc b/tests/rol/vector_adaptor_with_ghost_01.cc
new file mode 100644 (file)
index 0000000..50203d4
--- /dev/null
@@ -0,0 +1,159 @@
+// ---------------------------------------------------------------------
+//
+// 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 with MPI ghosted vectors using ROL::Vector's
+// checkVector() method.
+
+#include "../tests.h"
+
+#include <deal.II/lac/generic_linear_algebra.h>
+#include <deal.II/optimization/rol/vector_adaptor.h>
+
+using namespace dealii;
+
+// Vectors are prepared similar to 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 = 200;
+  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);
+  IndexSet local_relevant(global_size);
+  unsigned int ghost_indices [10] = {1, 2, 13, set-2, set-1, set, set+1, 2*set,
+                                     2*set+1, 2*set+3
+                                    };
+  local_relevant.add_indices (&ghost_indices[0], &ghost_indices[0]+10);
+
+  // --- Prepare vector.
+  v.reinit (local_owned, local_relevant, MPI_COMM_WORLD);
+}
+
+
+template <typename VectorType>
+void test ()
+{
+  VectorType a, b, c;
+  prepare_vector (a);
+  prepare_vector (b);
+  prepare_vector (c);
+
+  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;
+
+  for (auto iterator = c.begin(); iterator != c.end(); iterator++)
+    *iterator = static_cast<double>(Testing::rand())
+                /
+                RAND_MAX;
+
+  a.compress(VectorOperation::insert);
+  b.compress(VectorOperation::insert);
+  c.compress(VectorOperation::insert);
+
+  Teuchos::RCP<VectorType> a_rcp (new VectorType(a));
+  Teuchos::RCP<VectorType> b_rcp (new VectorType(b));
+  Teuchos::RCP<VectorType> c_rcp (new VectorType(c));
+
+  // --- Testing the constructor
+  Rol::VectorAdaptor<VectorType> a_rol (a_rcp);
+  Rol::VectorAdaptor<VectorType> b_rol (b_rcp);
+  Rol::VectorAdaptor<VectorType> c_rol (c_rcp);
+
+  Teuchos::RCP<std::ostream> out_stream;
+  Teuchos::oblackholestream bhs; // outputs nothing
+
+  if (dealii::Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0)
+    out_stream = Teuchos::rcp(&std::cout, false);
+  else
+    out_stream = Teuchos::rcp(&bhs, false);
+
+  a_rol.checkVector (b_rol, c_rol, true, *out_stream);
+}
+
+
+
+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<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_with_ghost_01.mpirun=10.output b/tests/rol/vector_adaptor_with_ghost_01.mpirun=10.output
new file mode 100644 (file)
index 0000000..c84a4b8
--- /dev/null
@@ -0,0 +1,18 @@
+
+********** Begin verification of linear algebra. *********************************************
+
+Commutativity of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Associativity of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Identity element of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Inverse elements of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Identity element of scalar multiplication. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Consistency of scalar multiplication with field multiplication. Consistency error: >>>>>>>>>>> 0
+Distributivity of scalar multiplication with respect to field addition. Consistency error: >>> 0
+Distributivity of scalar multiplication with respect to vector addition. Consistency error: >> 0
+Commutativity of dot (inner) product over the field of reals. Consistency error: >>>>>>>>>>>>> 0
+Additivity of dot (inner) product. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Consistency of scalar multiplication and norm. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Reflexivity. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+
+********** End verification of linear algebra. ***********************************************
+
diff --git a/tests/rol/vector_adaptor_with_ghost_01.mpirun=7.output b/tests/rol/vector_adaptor_with_ghost_01.mpirun=7.output
new file mode 100644 (file)
index 0000000..c84a4b8
--- /dev/null
@@ -0,0 +1,18 @@
+
+********** Begin verification of linear algebra. *********************************************
+
+Commutativity of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Associativity of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Identity element of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Inverse elements of addition. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Identity element of scalar multiplication. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Consistency of scalar multiplication with field multiplication. Consistency error: >>>>>>>>>>> 0
+Distributivity of scalar multiplication with respect to field addition. Consistency error: >>> 0
+Distributivity of scalar multiplication with respect to vector addition. Consistency error: >> 0
+Commutativity of dot (inner) product over the field of reals. Consistency error: >>>>>>>>>>>>> 0
+Additivity of dot (inner) product. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Consistency of scalar multiplication and norm. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+Reflexivity. Consistency error: >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 0
+
+********** End verification of linear algebra. ***********************************************
+

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.