Moved files from /reformulation to /ingredients/constraint_relaxation…

…_strategies, and from /solvers to /ingredients/subproblem_solvers

CMakeLists.txt

@@ -52,12 +52,12 @@ file(GLOB UNO_SOURCE_FILES
    uno/ingredients/inequality_handling_methods/*.cpp
    uno/ingredients/inequality_handling_methods/inequality_constrained_methods/*.cpp
    uno/ingredients/inequality_handling_methods/interior_point_methods/*.cpp
+   uno/ingredients/subproblem_solvers/*.cpp
+   uno/ingredients/subproblems/*.cpp
    uno/model/*.cpp
    uno/optimization/*.cpp
    uno/options/*.cpp
    uno/preprocessing/*.cpp
-   uno/reformulation/*.cpp
-   uno/solvers/*.cpp
    uno/tools/*.cpp
 )
 
@@ -115,14 +115,14 @@ else()
    endif()
 endif()
 if(HSL OR MA57)
-   list(APPEND UNO_SOURCE_FILES uno/solvers/MA57/MA57Solver.cpp)
+   list(APPEND UNO_SOURCE_FILES uno/ingredients/subproblem_solvers/MA57/MA57Solver.cpp)
    list(APPEND TESTS_UNO_SOURCE_FILES unotest/functional_tests/MA57SolverTests.cpp)
 
    find_package(BLAS REQUIRED)
    list(APPEND LIBRARIES ${BLAS_LIBRARIES})
 endif()
 if(HSL OR MA27)
-   list(APPEND UNO_SOURCE_FILES uno/solvers/MA27/MA27Solver.cpp)
+   list(APPEND UNO_SOURCE_FILES uno/ingredients/subproblem_solvers/MA27/MA27Solver.cpp)
    list(APPEND TESTS_UNO_SOURCE_FILES unotest/functional_tests/MA27SolverTests.cpp)
 
    find_package(BLAS REQUIRED)
@@ -142,7 +142,7 @@ find_library(BQPD bqpd)
 if(NOT BQPD)
    message(WARNING "Optional library BQPD was not found.")
 else()
-   list(APPEND UNO_SOURCE_FILES uno/solvers/BQPD/BQPDSolver.cpp)
+   list(APPEND UNO_SOURCE_FILES uno/ingredients/subproblem_solvers/BQPD/BQPDSolver.cpp)
    list(APPEND TESTS_UNO_SOURCE_FILES unotest/functional_tests/BQPDSolverTests.cpp)
    link_to_uno(bqpd ${BQPD})
 endif()
@@ -152,7 +152,7 @@ find_package(HIGHS)
 if(NOT HIGHS)
    message(WARNING "Optional library HiGHS was not found.")
 else()
-   list(APPEND UNO_SOURCE_FILES uno/solvers/HiGHS/HiGHSSolver.cpp)
+   list(APPEND UNO_SOURCE_FILES uno/ingredients/subproblem_solvers/HiGHS/HiGHSSolver.cpp)
    list(APPEND TESTS_UNO_SOURCE_FILES unotest/functional_tests/HiGHSSolverTests.cpp)
    link_to_uno(highs ${HIGHS})
    list(APPEND LIBRARIES highs::highs)
@@ -163,7 +163,7 @@ find_package(MUMPS)
 if(NOT MUMPS_LIBRARY)
    message(WARNING "Optional library MUMPS was not found.")
 else()
-   list(APPEND UNO_SOURCE_FILES uno/solvers/MUMPS/MUMPSSolver.cpp)
+   list(APPEND UNO_SOURCE_FILES uno/ingredients/subproblem_solvers/MUMPS/MUMPSSolver.cpp)
    list(APPEND TESTS_UNO_SOURCE_FILES unotest/functional_tests/MUMPSSolverTests.cpp)
    list(APPEND LIBRARIES ${MUMPS_LIBRARY} ${MUMPS_COMMON_LIBRARY} ${MUMPS_PORD_LIBRARY})
 

uno/Uno.cpp

@@ -8,13 +8,13 @@
 #include "ingredients/globalization_mechanisms/GlobalizationMechanismFactory.hpp"
 #include "ingredients/globalization_strategies/GlobalizationStrategyFactory.hpp"
 #include "ingredients/inequality_handling_methods/InequalityHandlingMethodFactory.hpp"
+#include "ingredients/subproblem_solvers/QPSolverFactory.hpp"
+#include "ingredients/subproblem_solvers/LPSolverFactory.hpp"
+#include "ingredients/subproblem_solvers/SymmetricIndefiniteLinearSolverFactory.hpp"
 #include "linear_algebra/Vector.hpp"
 #include "model/Model.hpp"
 #include "optimization/Iterate.hpp"
 #include "optimization/WarmstartInformation.hpp"
-#include "solvers/QPSolverFactory.hpp"
-#include "solvers/LPSolverFactory.hpp"
-#include "solvers/SymmetricIndefiniteLinearSolverFactory.hpp"
 #include "tools/Logger.hpp"
 #include "optimization/OptimizationStatus.hpp"
 #include "options/Options.hpp"

uno/ingredients/constraint_relaxation_strategies/ConstraintRelaxationStrategy.cpp

@@ -2,6 +2,7 @@
 // Licensed under the MIT license. See LICENSE file in the project directory for details.
 
 #include "ConstraintRelaxationStrategy.hpp"
+#include "OptimizationProblem.hpp"
 #include "ingredients/globalization_strategies/GlobalizationStrategy.hpp"
 #include "ingredients/globalization_strategies/GlobalizationStrategyFactory.hpp"
 #include "optimization/Direction.hpp"
@@ -11,7 +12,6 @@
 #include "model/Model.hpp"
 #include "optimization/Iterate.hpp"
 #include "optimization/Multipliers.hpp"
-#include "reformulation/OptimizationProblem.hpp"
 #include "symbolic/VectorView.hpp"
 #include "symbolic/Expression.hpp"
 #include "options/Options.hpp"

uno/reformulation/ElasticVariables.hpp → uno/ingredients/constraint_relaxation_strategies/ElasticVariables.hpp

@@ -30,4 +30,4 @@ namespace uno {
    };
 } // namespace
 
-#endif // UNO_ELASTICVARIABLES_H
+#endif // UNO_ELASTICVARIABLES_H

uno/ingredients/constraint_relaxation_strategies/FeasibilityRestoration.hpp

@@ -7,8 +7,8 @@
 #include <memory>
 #include "ConstraintRelaxationStrategy.hpp"
 #include "ingredients/globalization_strategies/ProgressMeasures.hpp"
-#include "reformulation/OptimalityProblem.hpp"
-#include "reformulation/l1RelaxedProblem.hpp"
+#include "OptimalityProblem.hpp"
+#include "l1RelaxedProblem.hpp"
 
 namespace uno {
    enum class Phase {FEASIBILITY_RESTORATION = 1, OPTIMALITY = 2};

uno/reformulation/OptimalityProblem.hpp → uno/ingredients/constraint_relaxation_strategies/OptimalityProblem.hpp

@@ -37,4 +37,4 @@ namespace uno {
    };
 } // namespace
 
-#endif // UNO_OPTIMALITYPROBLEM_H
+#endif // UNO_OPTIMALITYPROBLEM_H

uno/ingredients/constraint_relaxation_strategies/l1Relaxation.hpp

@@ -8,7 +8,7 @@
 #include "ConstraintRelaxationStrategy.hpp"
 #include "ingredients/globalization_strategies/ProgressMeasures.hpp"
 #include "optimization/Multipliers.hpp"
-#include "reformulation/l1RelaxedProblem.hpp"
+#include "l1RelaxedProblem.hpp"
 
 namespace uno {
    struct l1RelaxationParameters {

uno/ingredients/hessian_models/ConvexifiedHessian.cpp

@@ -4,12 +4,11 @@
 #include <stdexcept>
 #include "ConvexifiedHessian.hpp"
 #include "ingredients/hessian_models/UnstableRegularization.hpp"
-#include "reformulation/OptimizationProblem.hpp"
-#include "solvers/DirectSymmetricIndefiniteLinearSolver.hpp"
-#include "solvers/SymmetricIndefiniteLinearSolverFactory.hpp"
+#include "ingredients/constraint_relaxation_strategies/OptimizationProblem.hpp"
+#include "ingredients/subproblem_solvers/DirectSymmetricIndefiniteLinearSolver.hpp"
+#include "ingredients/subproblem_solvers/SymmetricIndefiniteLinearSolverFactory.hpp"
 #include "tools/Logger.hpp"
 #include "options/Options.hpp"
-#include "tools/Infinity.hpp"
 #include "tools/Statistics.hpp"
 
 namespace uno {

uno/ingredients/hessian_models/ExactHessian.cpp

@@ -2,7 +2,7 @@
 // Licensed under the MIT license. See LICENSE file in the project directory for details.
 
 #include "ExactHessian.hpp"
-#include "reformulation/OptimizationProblem.hpp"
+#include "ingredients/constraint_relaxation_strategies/OptimizationProblem.hpp"
 #include "options/Options.hpp"
 
 namespace uno {

uno/ingredients/hessian_models/HessianModelFactory.cpp

@@ -6,7 +6,7 @@
 #include "ConvexifiedHessian.hpp"
 #include "ExactHessian.hpp"
 #include "ZeroHessian.hpp"
-#include "solvers/DirectSymmetricIndefiniteLinearSolver.hpp"
+#include "ingredients/subproblem_solvers/DirectSymmetricIndefiniteLinearSolver.hpp"
 
 namespace uno {
    std::unique_ptr<HessianModel> HessianModelFactory::create(const std::string& hessian_model, size_t dimension, size_t maximum_number_nonzeros,

uno/ingredients/hessian_models/ZeroHessian.cpp

@@ -2,7 +2,7 @@
 // Licensed under the MIT license. See LICENSE file in the project directory for details.
 
 #include "ZeroHessian.hpp"
-#include "reformulation/OptimizationProblem.hpp"
+#include "ingredients/constraint_relaxation_strategies/OptimizationProblem.hpp"
 #include "options/Options.hpp"
 
 namespace uno {

uno/ingredients/inequality_handling_methods/InequalityHandlingMethodFactory.cpp

@@ -7,8 +7,8 @@
 #include "inequality_constrained_methods/QPSubproblem.hpp"
 #include "inequality_constrained_methods/LPSubproblem.hpp"
 #include "interior_point_methods/PrimalDualInteriorPointMethod.hpp"
-#include "solvers/QPSolverFactory.hpp"
-#include "solvers/SymmetricIndefiniteLinearSolverFactory.hpp"
+#include "ingredients/subproblem_solvers/QPSolverFactory.hpp"
+#include "ingredients/subproblem_solvers/SymmetricIndefiniteLinearSolverFactory.hpp"
 #include "options/Options.hpp"
 
 namespace uno {

uno/ingredients/inequality_handling_methods/inequality_constrained_methods/InequalityConstrainedMethod.cpp

@@ -5,7 +5,7 @@
 #include "optimization/Direction.hpp"
 #include "optimization/Iterate.hpp"
 #include "linear_algebra/Vector.hpp"
-#include "reformulation/l1RelaxedProblem.hpp"
+#include "ingredients/constraint_relaxation_strategies/l1RelaxedProblem.hpp"
 #include "options/Options.hpp"
 #include "symbolic/VectorView.hpp"
 

uno/ingredients/inequality_handling_methods/inequality_constrained_methods/LPSubproblem.cpp

@@ -4,9 +4,9 @@
 #include "LPSubproblem.hpp"
 #include "optimization/Direction.hpp"
 #include "optimization/WarmstartInformation.hpp"
-#include "reformulation/OptimizationProblem.hpp"
-#include "solvers/LPSolver.hpp"
-#include "solvers/LPSolverFactory.hpp"
+#include "ingredients/constraint_relaxation_strategies/OptimizationProblem.hpp"
+#include "ingredients/subproblem_solvers/LPSolver.hpp"
+#include "ingredients/subproblem_solvers/LPSolverFactory.hpp"
 #include "options/Options.hpp"
 
 namespace uno {

uno/ingredients/inequality_handling_methods/inequality_constrained_methods/QPSubproblem.cpp

@@ -3,15 +3,13 @@
 
 #include "QPSubproblem.hpp"
 #include "optimization/Direction.hpp"
-#include "ingredients/hessian_models/HessianModelFactory.hpp"
-#include "linear_algebra/SymmetricMatrix.hpp"
 #include "optimization/Iterate.hpp"
 #include "optimization/WarmstartInformation.hpp"
 #include "preprocessing/Preprocessing.hpp"
-#include "reformulation/OptimizationProblem.hpp"
-#include "solvers/DirectSymmetricIndefiniteLinearSolver.hpp"
-#include "solvers/QPSolver.hpp"
-#include "solvers/QPSolverFactory.hpp"
+#include "ingredients/constraint_relaxation_strategies/OptimizationProblem.hpp"
+#include "ingredients/subproblem_solvers/DirectSymmetricIndefiniteLinearSolver.hpp"
+#include "ingredients/subproblem_solvers/QPSolver.hpp"
+#include "ingredients/subproblem_solvers/QPSolverFactory.hpp"
 #include "options/Options.hpp"
 #include "tools/Statistics.hpp"
 

uno/ingredients/inequality_handling_methods/interior_point_methods/BarrierParameterUpdateStrategy.cpp

@@ -5,7 +5,7 @@
 #include <cmath>
 #include "BarrierParameterUpdateStrategy.hpp"
 #include "optimization/Iterate.hpp"
-#include "reformulation/OptimizationProblem.hpp"
+#include "ingredients/constraint_relaxation_strategies/OptimizationProblem.hpp"
 #include "symbolic/VectorExpression.hpp"
 #include "tools/Logger.hpp"
 #include "options/Options.hpp"

uno/ingredients/inequality_handling_methods/interior_point_methods/PrimalDualInteriorPointMethod.cpp

@@ -4,15 +4,14 @@
 #include <cmath>
 #include "PrimalDualInteriorPointMethod.hpp"
 #include "PrimalDualInteriorPointProblem.hpp"
+#include "ingredients/constraint_relaxation_strategies/l1RelaxedProblem.hpp"
+#include "ingredients/subproblem_solvers/DirectSymmetricIndefiniteLinearSolver.hpp"
+#include "ingredients/subproblem_solvers/SymmetricIndefiniteLinearSolverFactory.hpp"
+#include "linear_algebra/SparseStorageFactory.hpp"
 #include "optimization/Direction.hpp"
 #include "optimization/Iterate.hpp"
-#include "ingredients/hessian_models/HessianModelFactory.hpp"
-#include "linear_algebra/SparseStorageFactory.hpp"
-#include "solvers/DirectSymmetricIndefiniteLinearSolver.hpp"
-#include "solvers/SymmetricIndefiniteLinearSolverFactory.hpp"
 #include "optimization/WarmstartInformation.hpp"
 #include "preprocessing/Preprocessing.hpp"
-#include "reformulation/l1RelaxedProblem.hpp"
 #include "symbolic/VectorView.hpp"
 #include "tools/Infinity.hpp"
 
@@ -187,7 +186,7 @@ namespace uno {
 
       // check the inertia
       [[maybe_unused]] auto [number_pos_eigenvalues, number_neg_eigenvalues, number_zero_eigenvalues] = this->linear_solver->get_inertia();
-      assert(number_pos_eigenvalues == problem.number_variables && number_neg_eigenvalues == problem.number_constraints && number_zero_eigenvalues == 0);
+      assert(number_pos_eigenvalues == number_variables && number_neg_eigenvalues == number_constraints && number_zero_eigenvalues == 0);
 
       // rhs
       this->assemble_augmented_rhs(current_multipliers, number_variables, number_constraints);

uno/ingredients/inequality_handling_methods/interior_point_methods/PrimalDualInteriorPointProblem.hpp

@@ -4,7 +4,7 @@
 #ifndef UNO_PRIMALDUALINTERIORPOINTPROBLEM_H
 #define UNO_PRIMALDUALINTERIORPOINTPROBLEM_H
 
-#include "reformulation/OptimizationProblem.hpp"
+#include "ingredients/constraint_relaxation_strategies/OptimizationProblem.hpp"
 
 namespace uno {
    class PrimalDualInteriorPointProblem : public OptimizationProblem {

uno/solvers/BQPD/BQPDSolver.hpp → uno/ingredients/subproblem_solvers/BQPD/BQPDSolver.hpp

@@ -8,7 +8,7 @@
 #include <vector>
 #include "ingredients/subproblems/SubproblemStatus.hpp"
 #include "linear_algebra/Vector.hpp"
-#include "solvers/QPSolver.hpp"
+#include "ingredients/subproblem_solvers/QPSolver.hpp"
 
 namespace uno {
    // forward declaration

uno/ingredients/subproblem_solvers/BQPD/wdotd.f

@@ -0,0 +1,118 @@
+C Copyright (c) 2018-2024 Sven Leyffer
+C Licensed under the MIT license. See LICENSE file in the project directory for details.
+
+C     hristen this file wdotd.f
+
+      subroutine wdotd (n, x, ws, lws, result)
+
+c     ==========================================================
+c     Computes result = W.x where W is Hessian and x is a vector for AMPL
+c     Assumes v=0 on entry (OK, if called from gdotx, see QPsolve*.f)
+c     ==========================================================
+
+      implicit none
+
+c     ... declaration of passed parameters
+      integer          n, lws(0:*)
+      double precision x(n), result(n), ws(*)
+
+c     ... declaration of internal variables
+      integer i, j, k, footer_start
+
+c     inertia control for diagonal terms
+      double precision alpha
+      common /kktalphac/ alpha
+
+c     ========================  procedure body  =========================
+
+c     ... form result = W.x from sparse, upper triangular Hessian
+      footer_start = lws(0)
+      do i=1,n
+         do k=lws(footer_start + i - 1), lws(footer_start + i)-1
+            j = lws(k)
+            result(i) = result(i) + ws(k)*x(j)
+            if (j.ne.i) then
+c               off-diagonal term
+               result(j) = result(j) + ws(k)*x(i)
+            else
+c               diagonal term
+               result(i) = result(i) + alpha*x(j)
+            endif
+         enddo
+      enddo
+      return
+      end
+
+c     ******************************************************************
+
+      subroutine gdotx (n, x, ws, lws, result)
+
+      implicit none
+
+c     ... declaration of passed parameters
+      integer n, lws(*)
+      double precision    x(n), result(n), ws(*)
+
+c     ... declaration of internal variables
+      integer i
+
+c     ... storage map for hessian and scale_mode
+      integer         scale_mode, phe
+      common /scalec/ scale_mode, phe
+
+c     ========================  procedure body  =========================
+
+c     ... set result = 0
+      do i=1,n
+         result(i) = 0.D0
+      enddo
+
+c     ... allow for scaling of variables 
+      if ((scale_mode.eq.1).or.(scale_mode.eq.3)) then
+         do i=1,n
+            x(i) = x(i) * ws(i)
+         enddo
+      endif
+
+c     ... form v = W.d from sparse, upper triangular Hessian
+      call Wdotd (n, x, ws(phe+1), lws, result)
+
+c     ... allow for scaling of variables 
+      if ((scale_mode.eq.1).or.(scale_mode.eq.3)) then
+         do i=1,n
+            result(i) = result(i) * ws(i)
+            x(i) = x(i) / ws(i)
+         enddo
+      endif
+
+      return
+      end
+c     ******************************************************************
+      subroutine saipy2(s,a,la,i,y,n)
+      implicit double precision (a-h,o-z)
+      dimension a(*),la(0:*),y(*)
+c     ========================  procedure body  =========================
+c     saxpy with column i of A: y + s*A_{i, :}
+      if(s.eq.0.D0) return
+      j_column_start = la(0) + i
+      do j = la(j_column_start), la(j_column_start+1)-1
+         i_variable = la(j)
+         y(i_variable) = y(i_variable) + s*a(j)
+      enddo
+      return
+      end
+
+c     **************************** E N D *********************************
+      function daiscpr2(n,a,la,i,x,b)
+      implicit double precision (a-h,o-z)
+      dimension a(*),la(0:*),x(*)
+      DOUBLE PRECISION daiscpr2
+c     dot product of x and row i of A
+      daiscpr2 = dble(b)
+      j_column_start = la(0) + i
+      do j = la(j_column_start), la(j_column_start+1)-1
+         i_variable = la(j)
+         daiscpr2 = daiscpr2 + dble(x(i_variable))*dble(a(j))
+      enddo
+      return
+      end