]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Modernize Rol::VectorAdaptor to use ROL::Ptr. 17897/head
authorMarc Fehling <mafehling.git@gmail.com>
Thu, 28 Nov 2024 12:48:58 +0000 (13:48 +0100)
committerMarc Fehling <mafehling.git@gmail.com>
Thu, 28 Nov 2024 15:51:42 +0000 (16:51 +0100)
doc/news/changes/incompatibilities/20241128Fehling2 [new file with mode: 0644]
include/deal.II/optimization/rol/vector_adaptor.h
tests/rol/vector_adaptor_01.cc
tests/rol/vector_adaptor_02.cc
tests/rol/vector_adaptor_02.output
tests/rol/vector_adaptor_02.output.12.14.1 [deleted file]
tests/rol/vector_adaptor_02.output.12.4.2 [deleted file]
tests/rol/vector_adaptor_apply_binary_01.cc
tests/rol/vector_adaptor_no_ghost_01.cc
tests/rol/vector_adaptor_with_ghost_01.cc

diff --git a/doc/news/changes/incompatibilities/20241128Fehling2 b/doc/news/changes/incompatibilities/20241128Fehling2
new file mode 100644 (file)
index 0000000..74450bc
--- /dev/null
@@ -0,0 +1,6 @@
+Updated: The interface of Rol::VectorAdaptor class now uses `ROL::Ptr`.
+It is a shared pointer wrapper for either `Teuchos::RCP` or
+`std::shared_ptr`, and can be specified while configuring Trilinos.
+See `Trilinos/packages/rol/cmake/BuildOptions.cmake` for details.
+<br>
+(Marc Fehling, 2024/11/28)
index 2d50a3fbe145231665f8c3e83156cc06461886a2..d205a9901e63604fd6726936b7f4781bc0d4a6a2 100644 (file)
@@ -130,28 +130,26 @@ namespace Rol
 
   private:
     /**
-     * Teuchos smart reference counting pointer to the underlying vector of type
-     * <tt>VectorType</tt>.
+     * ROL pointer to the underlying vector of type <tt>VectorType</tt>.
      */
-    Teuchos::RCP<VectorType> vector_ptr;
+    ROL::Ptr<VectorType> vector_ptr;
 
   public:
     /**
      * Constructor.
      */
-    VectorAdaptor(const Teuchos::RCP<VectorType> &vector_ptr);
+    VectorAdaptor(const ROL::Ptr<VectorType> &vector_ptr);
 
     /**
-     * Return the Teuchos smart reference counting pointer to
-     * the wrapper vector, #vector_ptr.
+     * Return the ROL pointer to the wrapper vector, #vector_ptr.
      */
-    Teuchos::RCP<VectorType>
+    ROL::Ptr<VectorType>
     getVector();
 
     /**
-     * Return the Teuchos smart reference counting pointer to const vector.
+     * Return the ROL pointer to const vector.
      */
-    Teuchos::RCP<const VectorType>
+    ROL::Ptr<const VectorType>
     getVector() const;
 
     /**
@@ -211,15 +209,14 @@ namespace Rol
     /**
      * Return a clone of the wrapped vector.
      */
-    Teuchos::RCP<ROL::Vector<value_type>>
+    ROL::Ptr<ROL::Vector<value_type>>
     clone() const;
 
     /**
-     * Create and return a Teuchos smart reference counting pointer to the basis
-     * vector corresponding to the @p i ${}^{th}$ element of
-     * the wrapper vector.
+     * Create and return a ROL pointer to the basis vector corresponding to the
+     * @p i ${}^{th}$ element of the wrapper vector.
      */
-    Teuchos::RCP<ROL::Vector<value_type>>
+    ROL::Ptr<ROL::Vector<value_type>>
     basis(const int i) const;
 
     /**
@@ -257,14 +254,14 @@ namespace Rol
 
   template <typename VectorType>
   VectorAdaptor<VectorType>::VectorAdaptor(
-    const Teuchos::RCP<VectorType> &vector_ptr)
+    const ROL::Ptr<VectorType> &vector_ptr)
     : vector_ptr(vector_ptr)
   {}
 
 
 
   template <typename VectorType>
-  Teuchos::RCP<VectorType>
+  ROL::Ptr<VectorType>
   VectorAdaptor<VectorType>::getVector()
   {
     return vector_ptr;
@@ -273,7 +270,7 @@ namespace Rol
 
 
   template <typename VectorType>
-  Teuchos::RCP<const VectorType>
+  ROL::Ptr<const VectorType>
   VectorAdaptor<VectorType>::getVector() const
   {
     return vector_ptr;
@@ -286,7 +283,7 @@ namespace Rol
   VectorAdaptor<VectorType>::set(const ROL::Vector<value_type> &rol_vector)
   {
     const VectorAdaptor &vector_adaptor =
-      Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
+      dynamic_cast<const VectorAdaptor &>(rol_vector);
 
     (*vector_ptr) = *(vector_adaptor.getVector());
   }
@@ -301,7 +298,7 @@ namespace Rol
            ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
 
     const VectorAdaptor &vector_adaptor =
-      Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
+      dynamic_cast<const VectorAdaptor &>(rol_vector);
 
     *vector_ptr += *(vector_adaptor.getVector());
   }
@@ -317,7 +314,7 @@ namespace Rol
            ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
 
     const VectorAdaptor &vector_adaptor =
-      Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
+      dynamic_cast<const VectorAdaptor &>(rol_vector);
 
     vector_ptr->add(alpha, *(vector_adaptor.getVector()));
   }
@@ -354,7 +351,7 @@ namespace Rol
            ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
 
     const VectorAdaptor &vector_adaptor =
-      Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
+      dynamic_cast<const VectorAdaptor &>(rol_vector);
 
     return (*vector_ptr) * (*vector_adaptor.getVector());
   }
@@ -371,22 +368,21 @@ namespace Rol
 
 
   template <typename VectorType>
-  Teuchos::RCP<ROL::Vector<typename VectorType::value_type>>
+  ROL::Ptr<ROL::Vector<typename VectorType::value_type>>
   VectorAdaptor<VectorType>::clone() const
   {
-    Teuchos::RCP<VectorType> vec_ptr = Teuchos::rcp(new VectorType);
-    (*vec_ptr)                       = (*vector_ptr);
+    ROL::Ptr<VectorType> vec_ptr = ROL::makePtr<VectorType>(*vector_ptr);
 
-    return Teuchos::rcp(new VectorAdaptor(vec_ptr));
+    return ROL::makePtr<VectorAdaptor>(vec_ptr);
   }
 
 
 
   template <typename VectorType>
-  Teuchos::RCP<ROL::Vector<typename VectorType::value_type>>
+  ROL::Ptr<ROL::Vector<typename VectorType::value_type>>
   VectorAdaptor<VectorType>::basis(const int i) const
   {
-    Teuchos::RCP<VectorType> vec_ptr = Teuchos::rcp(new VectorType);
+    ROL::Ptr<VectorType> vec_ptr = ROL::makePtr<VectorType>();
 
     // Zero all the entries in dealii vector.
     vec_ptr->reinit(*vector_ptr, false);
@@ -403,9 +399,7 @@ namespace Rol
         vec_ptr->compress(VectorOperation::insert);
       }
 
-    Teuchos::RCP<VectorAdaptor> e = Teuchos::rcp(new VectorAdaptor(vec_ptr));
-
-    return e;
+    return ROL::makePtr<VectorAdaptor>(vec_ptr);
   }
 
 
@@ -444,7 +438,7 @@ namespace Rol
            ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
 
     const VectorAdaptor &vector_adaptor =
-      Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
+      dynamic_cast<const VectorAdaptor &>(rol_vector);
 
     const VectorType &given_rol_vector = *(vector_adaptor.getVector());
 
index af79567ab6f61c5d068c1cda3bc7ba03c89e73aa..9db7505bc2181134e8797fc3fc15b93ff4de2184 100644 (file)
@@ -25,7 +25,8 @@ template <typename VectorType>
 void
 test(const VectorType &given_vector)
 {
-  Teuchos::RCP<VectorType> given_vector_rcp(new VectorType(given_vector));
+  ROL::Ptr<VectorType> given_vector_rcp =
+    ROL::makePtr<VectorType>(given_vector);
 
   // --- Testing the constructor
   Rol::VectorAdaptor<VectorType> given_vector_rol(given_vector_rcp);
@@ -33,7 +34,7 @@ test(const VectorType &given_vector)
               ExcInternalError());
 
 
-  Teuchos::RCP<VectorType>       w_rcp = Teuchos::rcp(new VectorType);
+  ROL::Ptr<VectorType>           w_rcp = ROL::makePtr<VectorType>();
   Rol::VectorAdaptor<VectorType> w_rol(w_rcp);
 
   // --- Testing VectorAdaptor::set()
index 6f76c2e0cfe3b98f2e13fa61a5e42ace6a761399..4661792fa056d1406575f6db2bcca9b02c023d17 100644 (file)
@@ -24,7 +24,6 @@
 #include "ROL_LineSearchStep.hpp"
 #include "ROL_Objective.hpp"
 #include "ROL_StatusTest.hpp"
-#include "Teuchos_GlobalMPISession.hpp"
 
 // Use ROL to minimize the objective function, f(x,y) = x^2 + y^2.
 
@@ -36,21 +35,21 @@ 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)
+  ROL::Ptr<const VectorType>
+  get_rolptr_to_VectorType(const ROL::Vector<Real> &x)
   {
-    return (Teuchos::dyn_cast<const Xprim>(x)).getVector();
+    return (dynamic_cast<const Xprim &>(x)).getVector();
   }
 
-  Teuchos::RCP<dealii::Vector<Real>>
-  get_rcp_to_VectorType(ROL::Vector<Real> &x)
+  ROL::Ptr<dealii::Vector<Real>>
+  get_rolptr_to_VectorType(ROL::Vector<Real> &x)
   {
-    return (Teuchos::dyn_cast<Xprim>(x)).getVector();
+    return (dynamic_cast<Xprim &>(x)).getVector();
   }
 
 public:
   Real
-  value(const ROL::Vector<Real> &x, Real & /*tol*/)
+  value(const ROL::Vector<Real> &x, Real & /*tol*/) override
   {
     Assert(x.dimension() == 2, ExcInternalError());
 
@@ -58,10 +57,12 @@ public:
   }
 
   void
-  gradient(ROL::Vector<Real> &g, const ROL::Vector<Real> &x, Real & /*tol*/)
+  gradient(ROL::Vector<Real>       &g,
+           const ROL::Vector<Real> &x,
+           Real & /*tol*/) override
   {
-    Teuchos::RCP<const VectorType> xp = this->get_rcp_to_VectorType(x);
-    Teuchos::RCP<VectorType>       gp = this->get_rcp_to_VectorType(g);
+    ROL::Ptr<const VectorType> xp = this->get_rolptr_to_VectorType(x);
+    ROL::Ptr<VectorType>       gp = this->get_rolptr_to_VectorType(g);
 
     (*gp)[0] = 2. * (*xp)[0];
     (*gp)[1] = 2. * (*xp)[1];
@@ -75,8 +76,9 @@ test(const double x, const double y)
 
   QuadraticObjective<RealT> quad_objective;
 
-  Teuchos::RCP<std::ostream> outStream = Teuchos::rcp(&std::cout, false);
-  Teuchos::RCP<VectorType>   x_rcp     = Teuchos::rcp(new VectorType);
+  ROL::Ptr<std::ostream> outStream =
+    ROL::makePtrFromRef<std::ostream>(std::cout);
+  ROL::Ptr<VectorType> x_rcp = ROL::makePtr<VectorType>();
 
   x_rcp->reinit(2);
 
@@ -85,7 +87,7 @@ test(const double x, const double y)
 
   Rol::VectorAdaptor<VectorType> x_rol(x_rcp);
 
-  Teuchos::ParameterList parlist;
+  ROL::ParameterList parlist;
 
 #if DEAL_II_TRILINOS_VERSION_GTE(12, 18, 0)
   // Define algorithm in three intuitive and easy steps.
@@ -110,7 +112,7 @@ test(const double x, const double y)
   // Run Algorithm.
   algo.run(x_rol, quad_objective, true, *outStream);
 
-  Teuchos::RCP<const VectorType> xg = x_rol.getVector();
+  ROL::Ptr<const VectorType> xg = x_rol.getVector();
   std::cout << "The solution to minimization problem is: ";
   std::cout << (*xg)[0] << ' ' << (*xg)[1] << std::endl;
 }
index 75350d0a18aa4189d1ebb9f936c6b98d5a11de94..5c40863fde5149fea21cee8935a613d11d94d353 100644 (file)
@@ -3,19 +3,22 @@ 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         
+  1     0.000000e+00   0.000000e+00   1.019804e+01   3         2         2         0         
+Optimization Terminated with Status: Converged
 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         
+  1     0.000000e+00   0.000000e+00   1.414214e-01   3         2         2         0         
+Optimization Terminated with Status: Converged
 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         
+  1     0.000000e+00   0.000000e+00   1.095536e+01   3         2         2         0         
+Optimization Terminated with Status: Converged
 The solution to minimization problem is: 0 0
diff --git a/tests/rol/vector_adaptor_02.output.12.14.1 b/tests/rol/vector_adaptor_02.output.12.14.1
deleted file mode 100644 (file)
index 5c40863..0000000
+++ /dev/null
@@ -1,24 +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.040000e+02   2.039608e+01   
-  1     0.000000e+00   0.000000e+00   1.019804e+01   3         2         2         0         
-Optimization Terminated with Status: Converged
-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   3         2         2         0         
-Optimization Terminated with Status: Converged
-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   3         2         2         0         
-Optimization Terminated with Status: Converged
-The solution to minimization problem is: 0 0
diff --git a/tests/rol/vector_adaptor_02.output.12.4.2 b/tests/rol/vector_adaptor_02.output.12.4.2
deleted file mode 100644 (file)
index 2ae8dae..0000000
+++ /dev/null
@@ -1,21 +0,0 @@
-
-Quasi-Newton Method with Cubic Interpolation Linesearch satisfying Strong Wolfe Conditions
-Secant Type: Limited-Memory BFGS
-  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   3         2         2         0         
-The solution to minimization problem is: 0 0
-
-Quasi-Newton Method with Cubic Interpolation Linesearch satisfying Strong Wolfe Conditions
-Secant Type: Limited-Memory BFGS
-  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   3         2         2         0         
-The solution to minimization problem is: 0 0
-
-Quasi-Newton Method with Cubic Interpolation Linesearch satisfying Strong Wolfe Conditions
-Secant Type: Limited-Memory BFGS
-  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   3         2         2         0         
-The solution to minimization problem is: 0 0
index 13d030945bca4ccef0a46c10c6197d8ead73d275..768075887b7a35d00fe33f480fb2fd2332bbb4a4 100644 (file)
@@ -75,15 +75,15 @@ test()
   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::Ptr<VectorType> a_ptr = ROL::makePtr<VectorType>(a);
+  ROL::Ptr<VectorType> b_ptr = ROL::makePtr<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);
+  Rol::VectorAdaptor<VectorType> a_rol(a_ptr);
+  Rol::VectorAdaptor<VectorType> b_rol(b_ptr);
 
   a_rol.print(std::cout);
   b_rol.print(std::cout);
index 7eb23d67e090ef18c030f0e6fb5ac4c62c9c8b86..8724fd7dcecfc3d5f62a26e5b9aa159b43ce0823 100644 (file)
@@ -78,22 +78,22 @@ test()
   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));
+  ROL::Ptr<VectorType> a_ptr = ROL::makePtr<VectorType>(a);
+  ROL::Ptr<VectorType> b_ptr = ROL::makePtr<VectorType>(b);
+  ROL::Ptr<VectorType> c_ptr = ROL::makePtr<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);
+  Rol::VectorAdaptor<VectorType> a_rol(a_ptr);
+  Rol::VectorAdaptor<VectorType> b_rol(b_ptr);
+  Rol::VectorAdaptor<VectorType> c_rol(c_ptr);
 
-  Teuchos::RCP<std::ostream> out_stream;
-  Teuchos::oblackholestream  bhs; // outputs nothing
+  ROL::Ptr<std::ostream> out_stream;
+  ROL::nullstream        bhs; // outputs nothing
 
   if (dealii::Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0)
-    out_stream = Teuchos::rcp(&std::cout, false);
+    out_stream = ROL::makePtrFromRef<std::ostream>(std::cout);
   else
-    out_stream = Teuchos::rcp(&bhs, false);
+    out_stream = ROL::makePtrFromRef<std::ostream>(bhs);
 
   a_rol.checkVector(b_rol, c_rol, true, *out_stream);
 }
index 45c31b2db9d9574f8ae8b6850674a7ee823a90ed..a6796217c439cb298f4ceea8b1e83ff5a5fdda0c 100644 (file)
@@ -94,26 +94,26 @@ test()
   b.update_ghost_values();
   c.update_ghost_values();
 
-  Teuchos::RCP<VectorType> a_rcp(new VectorType(a));
-  Teuchos::RCP<VectorType> b_rcp(new VectorType(b));
-  Teuchos::RCP<VectorType> c_rcp(new VectorType(c));
+  ROL::Ptr<VectorType> a_ptr = ROL::makePtr<VectorType>(a);
+  ROL::Ptr<VectorType> b_ptr = ROL::makePtr<VectorType>(b);
+  ROL::Ptr<VectorType> c_ptr = ROL::makePtr<VectorType>(c);
 
-  a_rcp->update_ghost_values();
-  b_rcp->update_ghost_values();
-  c_rcp->update_ghost_values();
+  a_ptr->update_ghost_values();
+  b_ptr->update_ghost_values();
+  c_ptr->update_ghost_values();
 
   // --- Testing the constructor
-  Rol::VectorAdaptor<VectorType> a_rol(a_rcp);
-  Rol::VectorAdaptor<VectorType> b_rol(b_rcp);
-  Rol::VectorAdaptor<VectorType> c_rol(c_rcp);
+  Rol::VectorAdaptor<VectorType> a_rol(a_ptr);
+  Rol::VectorAdaptor<VectorType> b_rol(b_ptr);
+  Rol::VectorAdaptor<VectorType> c_rol(c_ptr);
 
-  Teuchos::RCP<std::ostream> out_stream;
-  Teuchos::oblackholestream  bhs; // outputs nothing
+  ROL::Ptr<std::ostream> out_stream;
+  ROL::nullstream        bhs; // outputs nothing
 
   if (dealii::Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0)
-    out_stream = Teuchos::rcp(&std::cout, false);
+    out_stream = ROL::makePtrFromRef<std::ostream>(std::cout);
   else
-    out_stream = Teuchos::rcp(&bhs, false);
+    out_stream = ROL::makePtrFromRef<std::ostream>(bhs);
 
   a_rol.checkVector(b_rol, c_rol, true, *out_stream);
 }

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.