]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add a vector adaptor to use Trilinos/ROL library; added cmake variable DEAL_II_TRILIN...
authorvishalkenchan <vishal.boddu@fau.de>
Fri, 10 Nov 2017 09:15:12 +0000 (10:15 +0100)
committervishalkenchan <vishal.boddu@fau.de>
Fri, 10 Nov 2017 12:27:59 +0000 (13:27 +0100)
cmake/configure/configure_2_trilinos.cmake
doc/news/changes/major/20171102VishalBoddu [new file with mode: 0644]
include/deal.II/base/config.h.in
include/deal.II/optimization/rol/vector_adaptor.h [new file with mode: 0644]
source/CMakeLists.txt
source/optimization/rol/CMakeLists.txt [new file with mode: 0644]

index deae1c811b000d8e06583601d572e3fdb4725fd9..d157125011e81e00ccbda85cceab1695bb1837bd 100644 (file)
@@ -54,6 +54,18 @@ MACRO(FEATURE_TRILINOS_FIND_EXTERNAL var)
       ENDIF()
     ENDFOREACH()
 
+    FOREACH(_optional_module
+         ROL
+      )
+      ITEM_MATCHES(_module_found ${_optional_module} ${Trilinos_PACKAGE_LIST})
+      IF(_module_found)
+        MESSAGE(STATUS "Found ${_optional_module}")
+        SET(DEAL_II_TRILINOS_WITH_${_optional_module} ON)
+      ELSE()
+        MESSAGE(STATUS "Module ${_optional_module} not found!")
+      ENDIF()
+    ENDFOREACH()
+
     IF((TRILINOS_VERSION_MAJOR EQUAL 11 AND
         NOT (TRILINOS_VERSION_MINOR LESS 14))
        OR
diff --git a/doc/news/changes/major/20171102VishalBoddu b/doc/news/changes/major/20171102VishalBoddu
new file mode 100644 (file)
index 0000000..c9c5d22
--- /dev/null
@@ -0,0 +1,7 @@
+New: The new namespace Rol contains an adaptor class that provides
+the implementation of the ROL::Vector interface. 
+The Trilinos package, Rapid Optimization Library (ROL), can solve unconstrained
+and constrained optimization problems, and optimization problems under
+uncertainty.
+<br>
+(Vishal Boddu 2017/11/02)
index fcf1065a62c94e89b8a551a20c511ee6fd3e0ab2..a19dffbc91d63bdb07c0f4e268292b8d3705928b 100644 (file)
 #cmakedefine DEAL_II_USE_MT_POSIX
 #cmakedefine DEAL_II_USE_MT_POSIX_NO_BARRIERS
 
+/* cmake/configure/configure_2_trilinos.cmake */
+#cmakedefine DEAL_II_TRILINOS_WITH_ROL
+
 /*
  * Depending on the use of threads, we will have to make some variables
  * volatile. We do this here in a very old-fashioned C-style, but still
diff --git a/include/deal.II/optimization/rol/vector_adaptor.h b/include/deal.II/optimization/rol/vector_adaptor.h
new file mode 100644 (file)
index 0000000..cc068a2
--- /dev/null
@@ -0,0 +1,489 @@
+//-----------------------------------------------------------
+//
+//    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.
+//
+//---------------------------------------------------------------
+
+#ifndef dealii_optimization_rol_vector_adaptor_h
+#define dealii_optimization_rol_vector_adaptor_h
+
+#include <deal.II/base/config.h>
+
+#ifdef DEAL_II_TRILINOS_WITH_ROL
+#include <ROL_Vector.hpp>
+
+#include <deal.II/base/exceptions.h>
+#include <deal.II/lac/vector.h>
+
+#include <type_traits>
+
+
+DEAL_II_NAMESPACE_OPEN
+
+
+using namespace dealii;
+
+
+/**
+ * A namespace that provides an interface to the
+ * <a href="https://trilinos.org/docs/dev/packages/rol/doc/html/index.html">
+ * Rapid Optimization Library</a> (ROL), a Trilinos package.
+ */
+namespace Rol
+{
+
+  /**
+   * An adaptor that provides the implementation of the ROL::Vector interface
+   * for vectors of type <tt>VectorType</tt>.
+   *
+   * This class supports vectors that satisfy the following requirements:
+   *
+   * The <tt>VectorType</tt> should contain the following types.
+   * ```
+   * VectorType::size_type;  // The type for size of the vector.
+   * VectorType::value_type; // The type for elements stored in the vector.
+   * VectorType::real_type;  // The type for real-valued numbers.
+   * ```
+   *
+   * However, ROL doesn't distinguish VectorAdaptor::value_type from
+   * VectorAdaptor::real_type. This is due to ROL's assumption that the
+   * VectorAdaptor::value_type itself is a type for real-valued numbers.
+   * Therefore, VectorAdaptor supports vectors whose real_type is
+   * convertible to value_type in the sense that
+   * <code>std::is_convertible<real_type, value_type>::value</code> yields
+   * <code>true</code>.
+   *
+   * The <tt>VectorType</tt> should contain the following methods.
+   * @code
+   * // Reinitialize the current vector using a given vector's
+   * // size (and the parallel distribution) without copying
+   * // the elements.
+   * VectorType::reinit(const VectorType &, ...);
+   *
+   * // Globally add a given vector to the current.
+   * VectorType::operator+=(const VectorType &);
+   *
+   * // Scale all elements by a given scalar.
+   * VectorType::operator*=(const VectorType::value_type &);
+   *
+   * // Perform dot product with a given vector.
+   * VectorType::operator*=(const VectorType &);
+   *
+   * // Scale all elements of the current vector and globally
+   * // add a given vector to it.
+   * VectorType::add(const VectorType::value_type, const VectorType &);
+   *
+   * // Copies the data of a given vector to the current.
+   * // Resize the current vector if necessary (MPI safe).
+   * VectorType::operation=(const VectorType &);
+   *
+   * // Return the global size of the current vector.
+   * VectorType::size();
+   *
+   * // Return L^2 norm of the current vector
+   * VectorType::l2_norm();
+   *
+   * // Iterator to the start of the (locally owned) element
+   * // of the current vector.
+   * VectorType::begin();
+   *
+   * // Iterator to the one past the last (locally owned)
+   * // element of the current vector.
+   * VectorType::end();
+   *
+   * // Compress the vector i.e., flush the buffers of the
+   * // vector object if it has any.
+   * VectorType::compress(VectorOperation::insert);
+   * @endcode
+   *
+   * @note The current implementation in ROL doesn't support vector sizes above
+   * the largest value of int type. Some of the vectors in deal.II (@ref Vector)
+   * may not satisfy the above requirements.
+   *
+   * @author Vishal Boddu, 2017
+   */
+  template<typename VectorType>
+  class VectorAdaptor : public ROL::Vector<typename VectorType::value_type>
+  {
+
+    /**
+     * A typedef for size type of <tt>VectorType</tt>.
+     */
+    using size_type = typename VectorType::size_type;
+
+    /**
+     * A typedef for element type stored in the <tt>VectorType</tt>.
+     */
+    using value_type = typename VectorType::value_type;
+
+    /**
+     * A typedef for real-valued numbers.
+     */
+    using real_type = typename VectorType::real_type;
+
+    static_assert (std::is_convertible<real_type, value_type>::value,
+                   "The real_type of the current VectorType is not "
+                   "convertible to the value_type.");
+
+  private:
+
+    /**
+     * Teuchos smart reference counting pointer to the underlying vector of type
+     * <tt>VectorType</tt>.
+     */
+    Teuchos::RCP<VectorType> vector_ptr;
+
+  public:
+
+    /**
+     * Constructor.
+     */
+    VectorAdaptor (const Teuchos::RCP<VectorType> &vector_ptr);
+
+    /**
+     * Returns the Teuchos smart reference counting pointer to
+     * the wrapper vector, #vector_ptr.
+     */
+    Teuchos::RCP<VectorType> getVector ();
+
+    /**
+     * Returns the Teuchos smart reference counting pointer to const vector.
+     */
+    Teuchos::RCP<const VectorType> getVector () const;
+
+    /**
+     * Return the dimension (global vector size) of the wrapped vector.
+     */
+    int dimension () const;
+
+    /**
+     * Set the wrapper vector to a given ROL::Vector @p rol_vector by
+     * overwriting its contents.
+     *
+     * If the current wrapper vector has ghost elements,
+     * then <code> VectorType::operator=(const VectorType&) </code> should still
+     * be allowed on it.
+     */
+    void set (const ROL::Vector<value_type> &rol_vector);
+
+    /**
+     * Perform addition.
+     */
+    void plus (const ROL::Vector<value_type> &rol_vector);
+
+    /**
+     * Scale the wrapper vector by @p alpha and add ROL::Vector @p rol_vector
+     * to it.
+     */
+    void axpy (const value_type               alpha,
+               const ROL::Vector<value_type> &rol_vector);
+
+    /**
+     * Scale the wrapper vector.
+     */
+    void scale (const value_type alpha);
+
+    /**
+     * Return the dot product with a given ROL::Vector @p rol_vector.
+     */
+    value_type dot (const ROL::Vector<value_type> &rol_vector) const;
+
+    /**
+     * Return the \f$ L^{2} \f$ norm of the wrapped vector.
+     *
+     * The returned type is of VectorAdaptor::value_type so as to maintain
+     * consistency with ROL::Vector<VectorAdaptor::value_type> and
+     * more importantly to not to create an overloaded version namely,
+     * <code> VectorAdaptor::real_type norm() const; </code>
+     * if real_type and value_type are not of the same type.
+     */
+    value_type norm() const;
+
+    /**
+     * Return a clone of the wrapped vector.
+     */
+    Teuchos::RCP<ROL::Vector<value_type>> clone() const;
+
+    /**
+     * Create and return a Teuchos smart reference counting pointer to the basis
+     * vector corresponding to the @p i \f${}^{th}\f$ element of
+     * the wrapper vector.
+     */
+    Teuchos::RCP<ROL::Vector<value_type>> basis (const int i) const;
+
+    /**
+     * Apply unary function @p f to all the elements of the wrapped vector.
+     */
+    void applyUnary (const ROL::Elementwise::UnaryFunction<value_type> &f);
+
+    /**
+     * 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);
+
+    /**
+     * Return the accumulated value on applying reduction operation @p r on
+     * all the elements of the wrapped vector.
+     */
+    value_type reduce (const ROL::Elementwise::ReductionOp<value_type> &r) const;
+
+    /**
+     * Print the wrapped vector to the output stream @p outStream.
+     */
+    void print (std::ostream &outStream) const;
+
+  };
+
+
+  /*------------------------------member definitions--------------------------*/
+#ifndef DOXYGEN
+
+
+  template<typename VectorType>
+  VectorAdaptor<VectorType>::
+  VectorAdaptor (const Teuchos::RCP<VectorType> &vector_ptr)
+    :
+    vector_ptr (vector_ptr)
+  {}
+
+
+
+  template<typename VectorType>
+  Teuchos::RCP<VectorType>
+  VectorAdaptor<VectorType>::getVector ()
+  {
+    return vector_ptr;
+  }
+
+
+
+  template<typename VectorType>
+  Teuchos::RCP<const VectorType>
+  VectorAdaptor<VectorType>::getVector () const
+  {
+    return vector_ptr;
+  }
+
+
+
+  template<typename VectorType>
+  void
+  VectorAdaptor<VectorType>::set (const ROL::Vector<value_type> &rol_vector)
+  {
+    const VectorAdaptor &vector_adaptor =
+      Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
+
+    (*vector_ptr) =  *(vector_adaptor.getVector());
+  }
+
+
+
+  template<typename VectorType>
+  void
+  VectorAdaptor<VectorType>::plus (const ROL::Vector<value_type> &rol_vector)
+  {
+    Assert (this->dimension() == rol_vector.dimension(),
+            ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
+
+    const VectorAdaptor &vector_adaptor =
+      Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
+
+    *vector_ptr += *(vector_adaptor.getVector());
+  }
+
+
+
+  template<typename VectorType>
+  void
+  VectorAdaptor<VectorType>::axpy (const value_type               alpha,
+                                   const ROL::Vector<value_type> &rol_vector)
+  {
+    Assert (this->dimension() == rol_vector.dimension(),
+            ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
+
+    const VectorAdaptor &vector_adaptor =
+      Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
+
+    vector_ptr->add (alpha, *(vector_adaptor.getVector()));
+  }
+
+
+
+  template<typename VectorType>
+  int
+  VectorAdaptor<VectorType>::dimension () const
+  {
+    Assert (vector_ptr->size() < std::numeric_limits<int>::max(),
+            ExcMessage("The size of the vector being used is greater than "
+                       "largest value of type int."));
+    return static_cast<int>(vector_ptr->size());
+  }
+
+
+
+  template<typename VectorType>
+  void
+  VectorAdaptor<VectorType>::scale (const value_type alpha)
+  {
+    (*vector_ptr) *= alpha;
+  }
+
+
+
+  template<typename VectorType>
+  typename VectorType::value_type
+  VectorAdaptor<VectorType>::
+  dot (const ROL::Vector<value_type> &rol_vector) const
+  {
+    Assert (this->dimension() == rol_vector.dimension(),
+            ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
+
+    const VectorAdaptor &vector_adaptor =
+      Teuchos::dyn_cast< const VectorAdaptor>(rol_vector);
+
+    return (*vector_ptr) * (*vector_adaptor.getVector());
+  }
+
+
+
+  template<typename VectorType>
+  typename VectorType::value_type
+  VectorAdaptor<VectorType>::norm() const
+  {
+    return vector_ptr->l2_norm();
+  }
+
+
+
+  template<typename VectorType>
+  Teuchos::RCP<ROL::Vector<typename VectorType::value_type> >
+  VectorAdaptor<VectorType>::clone() const
+  {
+    Teuchos::RCP< VectorType> vec_ptr = Teuchos::rcp (new VectorType);
+    (*vec_ptr) = (*vector_ptr);
+
+    return Teuchos::rcp (new VectorAdaptor (vec_ptr));
+  }
+
+
+
+  template<typename VectorType>
+  Teuchos::RCP<ROL::Vector<typename VectorType::value_type> >
+  VectorAdaptor<VectorType>::basis (const int i) const
+  {
+    Teuchos::RCP<VectorType> vec_ptr = Teuchos::rcp (new VectorType);
+
+    // Zero all the entries in dealii vector.
+    vec_ptr->reinit(*vector_ptr, false);
+
+    if (vector_ptr->locally_owned_elements().is_element(i))
+      vec_ptr->operator[](i) = 1.;
+
+    vec_ptr->compress(VectorOperation::insert);
+
+    Teuchos::RCP<VectorAdaptor> e = Teuchos::rcp (new VectorAdaptor(vec_ptr));
+
+    return e;
+  }
+
+
+
+  template<typename VectorType>
+  void
+  VectorAdaptor<VectorType>::
+  applyUnary (const ROL::Elementwise::UnaryFunction<value_type> &f)
+  {
+    const typename VectorType::iterator vend = vector_ptr->end();
+
+    for (typename VectorType::iterator
+         iterator  = vector_ptr->begin();
+         iterator != vend;
+         iterator++)
+      *iterator = f.apply(*iterator);
+
+    vector_ptr->compress (VectorOperation::insert);
+  }
+
+
+
+  template<typename VectorType>
+  void
+  VectorAdaptor<VectorType>::
+  applyBinary (const ROL::Elementwise::UnaryFunction<value_type> &f,
+               const ROL::Vector<value_type>                     &rol_vector)
+  {
+    Assert (this->dimension() == rol_vector.dimension(),
+            ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
+
+    const VectorAdaptor &vector_adaptor =
+      Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
+
+    const VectorType &given_rol_vector = *(vector_adaptor.getVector());
+
+    const typename VectorType::iterator
+    vend   = vector_ptr->end(),
+    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++)
+      *l_iterator = f.apply(*l_iterator, *r_iterator);
+
+    vector_ptr->compress (VectorOperation::insert);
+  }
+
+
+
+  template<typename VectorType>
+  typename VectorType::value_type
+  VectorAdaptor<VectorType>::
+  reduce (const ROL::Elementwise::ReductionOp<value_type> &r) const
+  {
+    typename VectorType::value_type result = r.initialValue();
+
+    const typename VectorType::iterator vend = vector_ptr->end();
+
+    for (typename VectorType::iterator
+         iterator  = vector_ptr->begin();
+         iterator != vend;
+         iterator++)
+      r.reduce(*iterator, result);
+    // Parallel reduction among processes is handled internally.
+
+    return result;
+  }
+
+
+
+  template<typename VectorType>
+  void
+  VectorAdaptor<VectorType>::print (std::ostream &outStream) const
+  {
+    vector_ptr->print(outStream);
+  }
+
+
+#endif // DOXYGEN
+
+
+} // namespace Rol
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+
+#endif // DEAL_II_TRILINOS_WITH_ROL
+
+#endif // dealii_optimization_rol_vector_adaptor_h
index d9c364a350532bce436455f59845acac6fdeb418..f789a8597f1cbb2ef3634fb8d711b7fc318f3d5d 100644 (file)
@@ -1,6 +1,6 @@
 ## ---------------------------------------------------------------------
 ##
-## Copyright (C) 2012 - 2016 by the deal.II authors
+## Copyright (C) 2012 - 2017 by the deal.II authors
 ##
 ## This file is part of the deal.II library.
 ##
@@ -50,6 +50,7 @@ ADD_SUBDIRECTORY(meshworker)
 ADD_SUBDIRECTORY(opencascade)
 ADD_SUBDIRECTORY(particles)
 ADD_SUBDIRECTORY(physics)
+ADD_SUBDIRECTORY(optimization/rol)
 ADD_SUBDIRECTORY(non_matching)
 ADD_SUBDIRECTORY(sundials)
 
diff --git a/source/optimization/rol/CMakeLists.txt b/source/optimization/rol/CMakeLists.txt
new file mode 100644 (file)
index 0000000..7cbeef9
--- /dev/null
@@ -0,0 +1,29 @@
+## ---------------------------------------------------------------------
+##
+## 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_DIRECTORIES(BEFORE ${CMAKE_CURRENT_BINARY_DIR})
+
+SET(_src
+  )
+
+SET(_inst
+  )
+
+FILE(GLOB _header
+  ${CMAKE_SOURCE_DIR}/include/deal.II/optimization/rol/*.h
+  )
+
+DEAL_II_ADD_LIBRARY(obj_rol OBJECT ${_src} ${_header} ${_inst})
+EXPAND_INSTANTIATIONS(obj_rol "${_inst}")

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.