From ad7d8e972625e26d8944dca27d3a8d57a6bf5987 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Thu, 28 Nov 2024 13:48:58 +0100 Subject: [PATCH] Modernize Rol::VectorAdaptor to use ROL::Ptr. --- .../incompatibilities/20241128Fehling2 | 6 ++ .../deal.II/optimization/rol/vector_adaptor.h | 56 +++++++++---------- tests/rol/vector_adaptor_01.cc | 5 +- tests/rol/vector_adaptor_02.cc | 32 ++++++----- tests/rol/vector_adaptor_02.output | 9 ++- tests/rol/vector_adaptor_02.output.12.14.1 | 24 -------- tests/rol/vector_adaptor_02.output.12.4.2 | 21 ------- tests/rol/vector_adaptor_apply_binary_01.cc | 8 +-- tests/rol/vector_adaptor_no_ghost_01.cc | 20 +++---- tests/rol/vector_adaptor_with_ghost_01.cc | 26 ++++----- 10 files changed, 84 insertions(+), 123 deletions(-) create mode 100644 doc/news/changes/incompatibilities/20241128Fehling2 delete mode 100644 tests/rol/vector_adaptor_02.output.12.14.1 delete mode 100644 tests/rol/vector_adaptor_02.output.12.4.2 diff --git a/doc/news/changes/incompatibilities/20241128Fehling2 b/doc/news/changes/incompatibilities/20241128Fehling2 new file mode 100644 index 0000000000..74450bc7d1 --- /dev/null +++ b/doc/news/changes/incompatibilities/20241128Fehling2 @@ -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. +
+(Marc Fehling, 2024/11/28) diff --git a/include/deal.II/optimization/rol/vector_adaptor.h b/include/deal.II/optimization/rol/vector_adaptor.h index 2d50a3fbe1..d205a9901e 100644 --- a/include/deal.II/optimization/rol/vector_adaptor.h +++ b/include/deal.II/optimization/rol/vector_adaptor.h @@ -130,28 +130,26 @@ namespace Rol private: /** - * Teuchos smart reference counting pointer to the underlying vector of type - * VectorType. + * ROL pointer to the underlying vector of type VectorType. */ - Teuchos::RCP vector_ptr; + ROL::Ptr vector_ptr; public: /** * Constructor. */ - VectorAdaptor(const Teuchos::RCP &vector_ptr); + VectorAdaptor(const ROL::Ptr &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 + ROL::Ptr getVector(); /** - * Return the Teuchos smart reference counting pointer to const vector. + * Return the ROL pointer to const vector. */ - Teuchos::RCP + ROL::Ptr getVector() const; /** @@ -211,15 +209,14 @@ namespace Rol /** * Return a clone of the wrapped vector. */ - Teuchos::RCP> + ROL::Ptr> 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::Ptr> basis(const int i) const; /** @@ -257,14 +254,14 @@ namespace Rol template VectorAdaptor::VectorAdaptor( - const Teuchos::RCP &vector_ptr) + const ROL::Ptr &vector_ptr) : vector_ptr(vector_ptr) {} template - Teuchos::RCP + ROL::Ptr VectorAdaptor::getVector() { return vector_ptr; @@ -273,7 +270,7 @@ namespace Rol template - Teuchos::RCP + ROL::Ptr VectorAdaptor::getVector() const { return vector_ptr; @@ -286,7 +283,7 @@ namespace Rol VectorAdaptor::set(const ROL::Vector &rol_vector) { const VectorAdaptor &vector_adaptor = - Teuchos::dyn_cast(rol_vector); + dynamic_cast(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(rol_vector); + dynamic_cast(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(rol_vector); + dynamic_cast(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(rol_vector); + dynamic_cast(rol_vector); return (*vector_ptr) * (*vector_adaptor.getVector()); } @@ -371,22 +368,21 @@ namespace Rol template - Teuchos::RCP> + ROL::Ptr> VectorAdaptor::clone() const { - Teuchos::RCP vec_ptr = Teuchos::rcp(new VectorType); - (*vec_ptr) = (*vector_ptr); + ROL::Ptr vec_ptr = ROL::makePtr(*vector_ptr); - return Teuchos::rcp(new VectorAdaptor(vec_ptr)); + return ROL::makePtr(vec_ptr); } template - Teuchos::RCP> + ROL::Ptr> VectorAdaptor::basis(const int i) const { - Teuchos::RCP vec_ptr = Teuchos::rcp(new VectorType); + ROL::Ptr vec_ptr = ROL::makePtr(); // 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 e = Teuchos::rcp(new VectorAdaptor(vec_ptr)); - - return e; + return ROL::makePtr(vec_ptr); } @@ -444,7 +438,7 @@ namespace Rol ExcDimensionMismatch(this->dimension(), rol_vector.dimension())); const VectorAdaptor &vector_adaptor = - Teuchos::dyn_cast(rol_vector); + dynamic_cast(rol_vector); const VectorType &given_rol_vector = *(vector_adaptor.getVector()); diff --git a/tests/rol/vector_adaptor_01.cc b/tests/rol/vector_adaptor_01.cc index af79567ab6..9db7505bc2 100644 --- a/tests/rol/vector_adaptor_01.cc +++ b/tests/rol/vector_adaptor_01.cc @@ -25,7 +25,8 @@ template void test(const VectorType &given_vector) { - Teuchos::RCP given_vector_rcp(new VectorType(given_vector)); + ROL::Ptr given_vector_rcp = + ROL::makePtr(given_vector); // --- Testing the constructor Rol::VectorAdaptor given_vector_rol(given_vector_rcp); @@ -33,7 +34,7 @@ test(const VectorType &given_vector) ExcInternalError()); - Teuchos::RCP w_rcp = Teuchos::rcp(new VectorType); + ROL::Ptr w_rcp = ROL::makePtr(); Rol::VectorAdaptor w_rol(w_rcp); // --- Testing VectorAdaptor::set() diff --git a/tests/rol/vector_adaptor_02.cc b/tests/rol/vector_adaptor_02.cc index 6f76c2e0cf..4661792fa0 100644 --- a/tests/rol/vector_adaptor_02.cc +++ b/tests/rol/vector_adaptor_02.cc @@ -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 QuadraticObjective : public ROL::Objective { private: - Teuchos::RCP - get_rcp_to_VectorType(const ROL::Vector &x) + ROL::Ptr + get_rolptr_to_VectorType(const ROL::Vector &x) { - return (Teuchos::dyn_cast(x)).getVector(); + return (dynamic_cast(x)).getVector(); } - Teuchos::RCP> - get_rcp_to_VectorType(ROL::Vector &x) + ROL::Ptr> + get_rolptr_to_VectorType(ROL::Vector &x) { - return (Teuchos::dyn_cast(x)).getVector(); + return (dynamic_cast(x)).getVector(); } public: Real - value(const ROL::Vector &x, Real & /*tol*/) + value(const ROL::Vector &x, Real & /*tol*/) override { Assert(x.dimension() == 2, ExcInternalError()); @@ -58,10 +57,12 @@ public: } void - gradient(ROL::Vector &g, const ROL::Vector &x, Real & /*tol*/) + gradient(ROL::Vector &g, + const ROL::Vector &x, + Real & /*tol*/) override { - Teuchos::RCP xp = this->get_rcp_to_VectorType(x); - Teuchos::RCP gp = this->get_rcp_to_VectorType(g); + ROL::Ptr xp = this->get_rolptr_to_VectorType(x); + ROL::Ptr 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 quad_objective; - Teuchos::RCP outStream = Teuchos::rcp(&std::cout, false); - Teuchos::RCP x_rcp = Teuchos::rcp(new VectorType); + ROL::Ptr outStream = + ROL::makePtrFromRef(std::cout); + ROL::Ptr x_rcp = ROL::makePtr(); x_rcp->reinit(2); @@ -85,7 +87,7 @@ test(const double x, const double y) Rol::VectorAdaptor 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 xg = x_rol.getVector(); + ROL::Ptr xg = x_rol.getVector(); std::cout << "The solution to minimization problem is: "; std::cout << (*xg)[0] << ' ' << (*xg)[1] << std::endl; } diff --git a/tests/rol/vector_adaptor_02.output b/tests/rol/vector_adaptor_02.output index 75350d0a18..5c40863fde 100644 --- a/tests/rol/vector_adaptor_02.output +++ b/tests/rol/vector_adaptor_02.output @@ -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 index 5c40863fde..0000000000 --- a/tests/rol/vector_adaptor_02.output.12.14.1 +++ /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 index 2ae8dae7a4..0000000000 --- a/tests/rol/vector_adaptor_02.output.12.4.2 +++ /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 diff --git a/tests/rol/vector_adaptor_apply_binary_01.cc b/tests/rol/vector_adaptor_apply_binary_01.cc index 13d030945b..768075887b 100644 --- a/tests/rol/vector_adaptor_apply_binary_01.cc +++ b/tests/rol/vector_adaptor_apply_binary_01.cc @@ -75,15 +75,15 @@ test() a.compress(VectorOperation::insert); b.compress(VectorOperation::insert); - Teuchos::RCP a_rcp(new VectorType(a)); - Teuchos::RCP b_rcp(new VectorType(b)); + ROL::Ptr a_ptr = ROL::makePtr(a); + ROL::Ptr b_ptr = ROL::makePtr(b); ROL::Elementwise::Multiply mult; ROL::Elementwise::Plus plus; // --- Testing the constructor - Rol::VectorAdaptor a_rol(a_rcp); - Rol::VectorAdaptor b_rol(b_rcp); + Rol::VectorAdaptor a_rol(a_ptr); + Rol::VectorAdaptor b_rol(b_ptr); a_rol.print(std::cout); b_rol.print(std::cout); diff --git a/tests/rol/vector_adaptor_no_ghost_01.cc b/tests/rol/vector_adaptor_no_ghost_01.cc index 7eb23d67e0..8724fd7dce 100644 --- a/tests/rol/vector_adaptor_no_ghost_01.cc +++ b/tests/rol/vector_adaptor_no_ghost_01.cc @@ -78,22 +78,22 @@ test() b.compress(VectorOperation::insert); c.compress(VectorOperation::insert); - Teuchos::RCP a_rcp(new VectorType(a)); - Teuchos::RCP b_rcp(new VectorType(b)); - Teuchos::RCP c_rcp(new VectorType(c)); + ROL::Ptr a_ptr = ROL::makePtr(a); + ROL::Ptr b_ptr = ROL::makePtr(b); + ROL::Ptr c_ptr = ROL::makePtr(c); // --- Testing the constructor - Rol::VectorAdaptor a_rol(a_rcp); - Rol::VectorAdaptor b_rol(b_rcp); - Rol::VectorAdaptor c_rol(c_rcp); + Rol::VectorAdaptor a_rol(a_ptr); + Rol::VectorAdaptor b_rol(b_ptr); + Rol::VectorAdaptor c_rol(c_ptr); - Teuchos::RCP out_stream; - Teuchos::oblackholestream bhs; // outputs nothing + ROL::Ptr 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::cout); else - out_stream = Teuchos::rcp(&bhs, false); + out_stream = ROL::makePtrFromRef(bhs); a_rol.checkVector(b_rol, c_rol, true, *out_stream); } diff --git a/tests/rol/vector_adaptor_with_ghost_01.cc b/tests/rol/vector_adaptor_with_ghost_01.cc index 45c31b2db9..a6796217c4 100644 --- a/tests/rol/vector_adaptor_with_ghost_01.cc +++ b/tests/rol/vector_adaptor_with_ghost_01.cc @@ -94,26 +94,26 @@ test() b.update_ghost_values(); c.update_ghost_values(); - Teuchos::RCP a_rcp(new VectorType(a)); - Teuchos::RCP b_rcp(new VectorType(b)); - Teuchos::RCP c_rcp(new VectorType(c)); + ROL::Ptr a_ptr = ROL::makePtr(a); + ROL::Ptr b_ptr = ROL::makePtr(b); + ROL::Ptr c_ptr = ROL::makePtr(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 a_rol(a_rcp); - Rol::VectorAdaptor b_rol(b_rcp); - Rol::VectorAdaptor c_rol(c_rcp); + Rol::VectorAdaptor a_rol(a_ptr); + Rol::VectorAdaptor b_rol(b_ptr); + Rol::VectorAdaptor c_rol(c_ptr); - Teuchos::RCP out_stream; - Teuchos::oblackholestream bhs; // outputs nothing + ROL::Ptr 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::cout); else - out_stream = Teuchos::rcp(&bhs, false); + out_stream = ROL::makePtrFromRef(bhs); a_rol.checkVector(b_rol, c_rol, true, *out_stream); } -- 2.39.5