Merge pull request #81 from cvanaret/feasibility_restoration

Improve the switch to optimality phase in FeasibilityRestoration

uno/ingredients/constraint_relaxation_strategies/ConstraintRelaxationStrategy.hpp

@@ -50,7 +50,7 @@ namespace uno {
 
       // trial iterate acceptance
       [[nodiscard]] virtual bool is_iterate_acceptable(Statistics& statistics, Iterate& current_iterate, Iterate& trial_iterate, const Direction& direction,
-            double step_length) = 0;
+            double step_length, WarmstartInformation& warmstart_information) = 0;
       [[nodiscard]] TerminationStatus check_termination(Iterate& iterate);
 
       // primal-dual residuals

uno/ingredients/constraint_relaxation_strategies/FeasibilityRestoration.cpp

@@ -118,12 +118,6 @@ namespace uno {
 
    void FeasibilityRestoration::solve_subproblem(Statistics& statistics, const OptimizationProblem& problem, Iterate& current_iterate,
          const Multipliers& current_multipliers, Direction& direction, WarmstartInformation& warmstart_information) {
-      // upon switching to the optimality phase, set a cold start in the subproblem solver
-      if (this->switching_to_optimality_phase) {
-         this->switching_to_optimality_phase = false;
-         warmstart_information.set_cold_start();
-      }
-
       direction.set_dimensions(problem.number_variables, problem.number_constraints);
       this->subproblem->solve(statistics, problem, current_iterate, current_multipliers, direction, warmstart_information);
       direction.norm = norm_inf(view(direction.primals, 0, this->model.number_variables));
@@ -138,7 +132,7 @@ namespace uno {
          direction.primals), this->residual_norm) <= this->linear_feasibility_tolerance);
    }
 
-   void FeasibilityRestoration::switch_to_optimality_phase(Iterate& current_iterate, Iterate& trial_iterate) {
+   void FeasibilityRestoration::switch_to_optimality_phase(Iterate& current_iterate, Iterate& trial_iterate, WarmstartInformation& warmstart_information) {
       DEBUG << "Switching from restoration to optimality phase\n";
       this->current_phase = Phase::OPTIMALITY;
       this->globalization_strategy->notify_switch_to_optimality(current_iterate.progress);
@@ -147,19 +141,20 @@ namespace uno {
       current_iterate.objective_multiplier = trial_iterate.objective_multiplier = 1.;
 
       this->subproblem->exit_feasibility_problem(this->optimality_problem, trial_iterate);
-      this->switching_to_optimality_phase = true;
+      // set a cold start in the subproblem solver
+      warmstart_information.set_cold_start();
    }
 
    bool FeasibilityRestoration::is_iterate_acceptable(Statistics& statistics, Iterate& current_iterate, Iterate& trial_iterate, const Direction& direction,
-         double step_length) {
+         double step_length, WarmstartInformation& warmstart_information) {
       // TODO pick right multipliers
       this->subproblem->postprocess_iterate(this->current_problem(), trial_iterate);
       this->compute_progress_measures(current_iterate, trial_iterate);
       trial_iterate.objective_multiplier = this->current_problem().get_objective_multiplier();
 
       // possibly go from restoration phase to optimality phase
       if (this->current_phase == Phase::FEASIBILITY_RESTORATION && this->can_switch_to_optimality_phase(current_iterate, trial_iterate, direction, step_length)) {
-         this->switch_to_optimality_phase(current_iterate, trial_iterate);
+         this->switch_to_optimality_phase(current_iterate, trial_iterate, warmstart_information);
       }
 
       bool accept_iterate = false;

uno/ingredients/constraint_relaxation_strategies/FeasibilityRestoration.hpp

@@ -29,7 +29,7 @@ namespace uno {
 
       // trial iterate acceptance
       [[nodiscard]] bool is_iterate_acceptable(Statistics& statistics, Iterate& current_iterate, Iterate& trial_iterate, const Direction& direction,
-            double step_length) override;
+            double step_length, WarmstartInformation& warmstart_information) override;
 
       // primal-dual residuals
       void compute_primal_dual_residuals(Iterate& iterate) override;
@@ -41,7 +41,6 @@ namespace uno {
       Phase current_phase{Phase::OPTIMALITY};
       const double linear_feasibility_tolerance;
       const bool switch_to_optimality_requires_linearized_feasibility;
-      bool switching_to_optimality_phase{false};
       ProgressMeasures reference_optimality_progress{};
       Vector<double> reference_optimality_primals{};
 
@@ -51,7 +50,7 @@ namespace uno {
       [[nodiscard]] const OptimizationProblem& current_problem() const;
       void solve_subproblem(Statistics& statistics, const OptimizationProblem& problem, Iterate& current_iterate, const Multipliers& current_multipliers,
             Direction& direction, WarmstartInformation& warmstart_information);
-      void switch_to_optimality_phase(Iterate& current_iterate, Iterate& trial_iterate);
+      void switch_to_optimality_phase(Iterate& current_iterate, Iterate& trial_iterate, WarmstartInformation& warmstart_information);
 
       void evaluate_progress_measures(Iterate& iterate) const override;
       [[nodiscard]] ProgressMeasures compute_predicted_reduction_models(Iterate& current_iterate, const Direction& direction, double step_length);

uno/ingredients/constraint_relaxation_strategies/l1Relaxation.cpp

@@ -232,7 +232,7 @@ namespace uno {
    }
 
    bool l1Relaxation::is_iterate_acceptable(Statistics& statistics, Iterate& current_iterate, Iterate& trial_iterate, const Direction& direction,
-         double step_length) {
+         double step_length, WarmstartInformation& /*warmstart_information*/) {
       this->subproblem->postprocess_iterate(this->l1_relaxed_problem, trial_iterate);
       this->compute_progress_measures(current_iterate, trial_iterate);
       trial_iterate.objective_multiplier = this->l1_relaxed_problem.get_objective_multiplier();

uno/ingredients/constraint_relaxation_strategies/l1Relaxation.hpp

@@ -36,7 +36,7 @@ namespace uno {
 
       // trial iterate acceptance
       [[nodiscard]] bool is_iterate_acceptable(Statistics& statistics, Iterate& current_iterate, Iterate& trial_iterate, const Direction& direction,
-            double step_length) override;
+            double step_length, WarmstartInformation& warmstart_information) override;
 
       // primal-dual residuals
       void compute_primal_dual_residuals(Iterate& iterate) override;

uno/ingredients/globalization_mechanisms/BacktrackingLineSearch.cpp

@@ -59,7 +59,8 @@ namespace uno {
                   // scale or not the constraint dual direction with the LS step length
                   this->scale_duals_with_step_length ? step_length : 1.);
 
-            is_acceptable = this->constraint_relaxation_strategy.is_iterate_acceptable(statistics, current_iterate, trial_iterate, this->direction, step_length);
+            is_acceptable = this->constraint_relaxation_strategy.is_iterate_acceptable(statistics, current_iterate, trial_iterate, this->direction,
+                  step_length, warmstart_information);
             this->set_statistics(statistics, trial_iterate, this->direction, step_length, number_iterations);
          }
          catch (const EvaluationError& e) {

uno/ingredients/globalization_mechanisms/TrustRegionStrategy.cpp

@@ -78,7 +78,7 @@ namespace uno {
                GlobalizationMechanism::assemble_trial_iterate(model, current_iterate, trial_iterate, this->direction, 1., 1.);
                this->reset_active_trust_region_multipliers(model, this->direction, trial_iterate);
 
-               is_acceptable = this->is_iterate_acceptable(statistics, current_iterate, trial_iterate, this->direction);
+               is_acceptable = this->is_iterate_acceptable(statistics, current_iterate, trial_iterate, this->direction, warmstart_information);
                if (is_acceptable) {
                   this->constraint_relaxation_strategy.set_dual_residuals_statistics(statistics, trial_iterate);
                   this->reset_radius();
@@ -124,9 +124,10 @@ namespace uno {
 
    // the trial iterate is accepted by the constraint relaxation strategy or if the step is small and we cannot switch to solving the feasibility problem
    bool TrustRegionStrategy::is_iterate_acceptable(Statistics& statistics, Iterate& current_iterate, Iterate& trial_iterate,
-         const Direction& direction) {
+         const Direction& direction, WarmstartInformation& warmstart_information) {
       // direction.primal_dual_step_length is usually 1, can be lower if reduced by fraction-to-boundary rule
-      bool accept_iterate = this->constraint_relaxation_strategy.is_iterate_acceptable(statistics, current_iterate, trial_iterate, direction, 1.);
+      bool accept_iterate = this->constraint_relaxation_strategy.is_iterate_acceptable(statistics, current_iterate, trial_iterate, direction, 1.,
+            warmstart_information);
       this->set_statistics(statistics, trial_iterate, direction);
       if (accept_iterate) {
          trial_iterate.status = this->constraint_relaxation_strategy.check_termination(trial_iterate);

uno/ingredients/globalization_mechanisms/TrustRegionStrategy.hpp

@@ -24,7 +24,8 @@ namespace uno {
       const double radius_reset_threshold;
       const double tolerance;
 
-      bool is_iterate_acceptable(Statistics& statistics, Iterate& current_iterate, Iterate& trial_iterate, const Direction& direction);
+      bool is_iterate_acceptable(Statistics& statistics, Iterate& current_iterate, Iterate& trial_iterate, const Direction& direction,
+            WarmstartInformation& warmstart_information);
       void possibly_increase_radius(double step_norm);
       void decrease_radius(double step_norm);
       void decrease_radius();