Skip to content

Commit

Permalink
Merge branch 'latest' into latest-merge
Browse files Browse the repository at this point in the history
  • Loading branch information
@galabovaa
galabovaa committed May 26, 2023
2 parents 1062cfb + a0c6768 commit 975dfc3
Show file tree
Hide file tree
Showing 11 changed files with 104 additions and 68 deletions.
12 changes: 6 additions & 6 deletions extern/filereaderlp/reader.cpp
View file
Expand Up @@ -828,7 +828,7 @@ void Reader::splittokens() {
// Current section is non-trivial, so mark its end, using the
// value of currentsection to indicate that there is no open
// section
assert(debug_open_section);
lpassert(debug_open_section);
sectiontokens[currentsection].second = it;
debug_open_section = false;
currentsection = LpSectionKeyword::NONE;
Expand All @@ -846,12 +846,12 @@ void Reader::splittokens() {
// the current section
if (currentsection != LpSectionKeyword::NONE &&
currentsection != next->keyword) {
assert(debug_open_section);
lpassert(debug_open_section);
sectiontokens[currentsection].second = it;
debug_open_section = false;
}
currentsection = LpSectionKeyword::NONE;
assert(!debug_open_section);
lpassert(!debug_open_section);
continue;
}
// Next section is non-empty
Expand All @@ -862,16 +862,16 @@ void Reader::splittokens() {
lpassert(sectiontokens.count(currentsection) == 0);
// Remember the beginning of the new section: its the token
// following the current one
assert(!debug_open_section);
lpassert(!debug_open_section);
sectiontokens[currentsection].first = next;
debug_open_section = true;
}
// Always ends with either an open section or a section type of
// LpSectionKeyword::NONE
assert(debug_open_section != (currentsection == LpSectionKeyword::NONE));
lpassert(debug_open_section != (currentsection == LpSectionKeyword::NONE));
}
// Check that the last section has been closed
assert(currentsection == LpSectionKeyword::NONE);
lpassert(currentsection == LpSectionKeyword::NONE);
}

void Reader::processtokens() {
Expand Down
41 changes: 36 additions & 5 deletions src/lp_data/Highs.cpp
View file
Expand Up @@ -3150,17 +3150,18 @@ HighsStatus Highs::callSolveQp() {
return_status = interpretCallStatus(options_.log_options, call_status,
return_status, "QpSolver");
if (return_status == HighsStatus::kError) return return_status;
model_status_ = runtime.status == ProblemStatus::OPTIMAL
model_status_ = runtime.status == QpModelStatus::OPTIMAL
? HighsModelStatus::kOptimal
: runtime.status == ProblemStatus::UNBOUNDED
: runtime.status == QpModelStatus::UNBOUNDED
? HighsModelStatus::kUnbounded
: runtime.status == ProblemStatus::INFEASIBLE
: runtime.status == QpModelStatus::INFEASIBLE
? HighsModelStatus::kInfeasible
: runtime.status == ProblemStatus::ITERATIONLIMIT
: runtime.status == QpModelStatus::ITERATIONLIMIT
? HighsModelStatus::kIterationLimit
: runtime.status == ProblemStatus::TIMELIMIT
: runtime.status == QpModelStatus::TIMELIMIT
? HighsModelStatus::kTimeLimit
: HighsModelStatus::kNotset;
// extract variable values
solution_.col_value.resize(lp.num_col_);
solution_.col_dual.resize(lp.num_col_);
const double objective_multiplier = lp.sense_ == ObjSense::kMinimize ? 1 : -1;
Expand All @@ -3169,6 +3170,7 @@ HighsStatus Highs::callSolveQp() {
solution_.col_dual[iCol] =
objective_multiplier * runtime.dualvar.value[iCol];
}
// extract constraint activity
solution_.row_value.resize(lp.num_row_);
solution_.row_dual.resize(lp.num_row_);
// Negate the vector and Hessian
Expand All @@ -3179,6 +3181,35 @@ HighsStatus Highs::callSolveQp() {
}
solution_.value_valid = true;
solution_.dual_valid = true;

// extract basis status
basis_.col_status.resize(lp.num_col_);
basis_.row_status.resize(lp.num_row_);

for (HighsInt i = 0; i < lp.num_col_; i++) {
if (runtime.status_var[i] == BasisStatus::ActiveAtLower) {
basis_.col_status[i] = HighsBasisStatus::kLower;
} else if (runtime.status_var[i] == BasisStatus::ActiveAtUpper) {
basis_.col_status[i] = HighsBasisStatus::kUpper;
} else if (runtime.status_var[i] == BasisStatus::InactiveInBasis) {
basis_.col_status[i] = HighsBasisStatus::kNonbasic;
} else {
basis_.col_status[i] = HighsBasisStatus::kBasic;
}
}

for (HighsInt i = 0; i < lp.num_row_; i++) {
if (runtime.status_con[i] == BasisStatus::ActiveAtLower) {
basis_.row_status[i] = HighsBasisStatus::kLower;
} else if (runtime.status_con[i] == BasisStatus::ActiveAtUpper) {
basis_.row_status[i] = HighsBasisStatus::kUpper;
} else if (runtime.status_con[i] == BasisStatus::InactiveInBasis) {
basis_.row_status[i] = HighsBasisStatus::kNonbasic;
} else {
basis_.row_status[i] = HighsBasisStatus::kBasic;
}
}

// Get the objective and any KKT failures
info_.objective_function_value = model_.objectiveValue(solution_.col_value);
getKktFailures(options_, model_, solution_, basis_, info_);
Expand Down
4 changes: 2 additions & 2 deletions src/presolve/HighsPostsolveStack.cpp
View file
Expand Up @@ -568,7 +568,7 @@ void HighsPostsolveStack::DuplicateRow::undo(const HighsOptions& options,
void HighsPostsolveStack::DuplicateColumn::undo(const HighsOptions& options,
HighsSolution& solution,
HighsBasis& basis) const {
const bool debug_report = true;
const bool debug_report = false;
const double mergeVal = solution.col_value[col];

auto okResidual = [&](const double x, const double y) {
Expand Down Expand Up @@ -986,7 +986,7 @@ void HighsPostsolveStack::DuplicateColumn::undoFix(
options.primal_feasibility_tolerance;
std::vector<double>& col_value = solution.col_value;
const bool allow_assert = false;
const bool debug_report = true;
const bool debug_report = false;
//=============================================================================================

auto isInteger = [&](const double v) {
Expand Down
21 changes: 14 additions & 7 deletions src/qpsolver/basis.cpp
View file
Expand Up @@ -4,17 +4,23 @@
#include <memory>

Basis::Basis(Runtime& rt, std::vector<HighsInt> active,
std::vector<BasisStatus> lower, std::vector<HighsInt> inactive)
std::vector<BasisStatus> status, std::vector<HighsInt> inactive)
: runtime(rt),
buffer_column_aq(rt.instance.num_var),
buffer_row_ep(rt.instance.num_var) {
buffer_vec2hvec.setup(rt.instance.num_var);

for (HighsInt i=0; i<runtime.instance.num_var + runtime.instance.num_con; i++) {
basisstatus[i] = BasisStatus::Inactive;
}

for (size_t i = 0; i < active.size(); i++) {
activeconstraintidx.push_back(active[i]);
basisstatus[activeconstraintidx[i]] = lower[i];
basisstatus[activeconstraintidx[i]] = status[i];
}
for (HighsInt i : inactive) {
nonactiveconstraintsidx.push_back(i);
for (size_t i=0; i<inactive.size(); i++) {
nonactiveconstraintsidx.push_back(inactive[i]);
basisstatus[nonactiveconstraintsidx[i]] = BasisStatus::InactiveInBasis;
}

Atran = rt.instance.A.t();
Expand Down Expand Up @@ -99,16 +105,17 @@ void Basis::report() {
void Basis::deactivate(HighsInt conid) {
// printf("deact %" HIGHSINT_FORMAT "\n", conid);
assert(contains(activeconstraintidx, conid));
basisstatus.erase(conid);
basisstatus[conid] = BasisStatus::InactiveInBasis;
remove(activeconstraintidx, conid);
nonactiveconstraintsidx.push_back(conid);
}

QpSolverStatus Basis::activate(const Settings& settings, HighsInt conid, BasisStatus atlower,
QpSolverStatus Basis::activate(const Settings& settings, HighsInt conid, BasisStatus newstatus,
HighsInt nonactivetoremove, Pricing* pricing) {
// printf("activ %" HIGHSINT_FORMAT "\n", conid);
if (!contains(activeconstraintidx, (HighsInt)conid)) {
basisstatus[conid] = atlower;
basisstatus[nonactivetoremove] = BasisStatus::Inactive;
basisstatus[conid] = newstatus;
activeconstraintidx.push_back(conid);
} else {
printf("Degeneracy? constraint %" HIGHSINT_FORMAT " already in basis\n",
Expand Down
8 changes: 0 additions & 8 deletions src/qpsolver/basis.hpp
View file
Expand Up @@ -14,14 +14,6 @@
#include "util/HVector.h"
#include "util/HVectorBase.h"

enum class BasisStatus {
Default,
ActiveAtLower = 1,
ActiveAtUpper,
ActiveAtZero,
Inactive
};

class Basis {
HVector buffer_vec2hvec;

Expand Down
4 changes: 2 additions & 2 deletions src/qpsolver/feasibility_highs.hpp
View file
Expand Up @@ -70,15 +70,15 @@ static void computestartingpoint_highs(Runtime& runtime, CrashSolution& result)

HighsStatus status = highs.run();
if (status != HighsStatus::kOk) {
runtime.status = ProblemStatus::ERROR;
runtime.status = QpModelStatus::ERROR;
return;
}

runtime.statistics.phase1_iterations = highs.getInfo().simplex_iteration_count;

HighsModelStatus phase1stat = highs.getModelStatus();
if (phase1stat == HighsModelStatus::kInfeasible) {
runtime.status = ProblemStatus::INFEASIBLE;
runtime.status = QpModelStatus::INFEASIBLE;
return;
}

Expand Down
2 changes: 1 addition & 1 deletion src/qpsolver/feasibility_quass.hpp
View file
Expand Up @@ -57,7 +57,7 @@ void computestartingpoint_quass(Runtime& runtime, CrashSolution& result) {

// HighsModelStatus phase1stat = highs.getModelStatus();
// if (phase1stat == HighsModelStatus::kInfeasible) {
// runtime.status = ProblemStatus::INFEASIBLE;
// runtime.status = QpModelStatus::INFEASIBLE;
// return;
// }

Expand Down
18 changes: 18 additions & 0 deletions src/qpsolver/qpconst.hpp
View file
Expand Up @@ -4,4 +4,22 @@

enum class QpSolverStatus { OK, NOTPOSITIVDEFINITE, DEGENERATE };

enum class QpModelStatus {
INDETERMINED,
OPTIMAL,
UNBOUNDED,
INFEASIBLE,
ITERATIONLIMIT,
TIMELIMIT,
ERROR
};

enum class BasisStatus {
Inactive,
ActiveAtLower = 1,
ActiveAtUpper,
InactiveInBasis
};


#endif
25 changes: 17 additions & 8 deletions src/qpsolver/quass.cpp
View file
Expand Up @@ -29,7 +29,7 @@ void Quass::solve() {
runtime.instance = runtime.perturbed;
CrashSolution crash(runtime.instance.num_var, runtime.instance.num_con);
computestartingpoint(runtime, crash);
if (runtime.status != ProblemStatus::INDETERMINED) {
if (runtime.status != QpModelStatus::INDETERMINED) {
return;
}
Basis basis(runtime, crash.active, crash.rowstatus, crash.inactive);
Expand Down Expand Up @@ -246,13 +246,13 @@ void Quass::solve(const Vector& x0, const Vector& ra, Basis& b0) {
while (true) {
// check iteration limit
if (runtime.statistics.num_iterations >= runtime.settings.iterationlimit) {
runtime.status = ProblemStatus::ITERATIONLIMIT;
runtime.status = QpModelStatus::ITERATIONLIMIT;
break;
}

// check time limit
if (runtime.timer.readRunHighsClock() >= runtime.settings.timelimit) {
runtime.status = ProblemStatus::TIMELIMIT;
runtime.status = QpModelStatus::TIMELIMIT;
break;
}

Expand All @@ -272,7 +272,7 @@ void Quass::solve(const Vector& x0, const Vector& ra, Basis& b0) {
if (atfsep) {
HighsInt minidx = pricing->price(runtime.primal, gradient.getGradient());
if (minidx == -1) {
runtime.status = ProblemStatus::OPTIMAL;
runtime.status = QpModelStatus::OPTIMAL;
break;
}

Expand All @@ -295,7 +295,7 @@ void Quass::solve(const Vector& x0, const Vector& ra, Basis& b0) {
if (!zero_curvature_direction) {
status = factor.expand(buffer_yp, buffer_gyp, buffer_l, buffer_m);
if (status != QpSolverStatus::OK) {
runtime.status = ProblemStatus::INDETERMINED;
runtime.status = QpModelStatus::INDETERMINED;
return;
}
}
Expand All @@ -317,7 +317,7 @@ void Quass::solve(const Vector& x0, const Vector& ra, Basis& b0) {
status = reduce(runtime, basis, stepres.limitingconstraint, buffer_d,
maxabsd, constrainttodrop);
if (status != QpSolverStatus::OK) {
runtime.status = ProblemStatus::INDETERMINED;
runtime.status = QpModelStatus::INDETERMINED;
return;
}
if (!zero_curvature_direction) {
Expand All @@ -334,7 +334,7 @@ void Quass::solve(const Vector& x0, const Vector& ra, Basis& b0) {
: BasisStatus::ActiveAtUpper,
constrainttodrop, pricing.get());
if (status != QpSolverStatus::OK) {
runtime.status = ProblemStatus::INDETERMINED;
runtime.status = QpModelStatus::INDETERMINED;
return;
}
if (basis.getnumactive() != runtime.instance.num_var) {
Expand All @@ -344,7 +344,7 @@ void Quass::solve(const Vector& x0, const Vector& ra, Basis& b0) {
if (stepres.alpha ==
std::numeric_limits<double>::infinity()) {
// unbounded
runtime.status = ProblemStatus::UNBOUNDED;
runtime.status = QpModelStatus::UNBOUNDED;
return;
}
atfsep = false;
Expand Down Expand Up @@ -377,6 +377,15 @@ void Quass::solve(const Vector& x0, const Vector& ra, Basis& b0) {
}
}

// extract basis status
for (HighsInt i=0; i<runtime.instance.num_var; i++) {
runtime.status_var[i] = basis.getstatus(i);
}

for (HighsInt i=0; i<runtime.instance.num_con; i++) {
runtime.status_con[i] = basis.getstatus(runtime.instance.num_var + i);
}

if (basis.getnumactive() == runtime.instance.num_var) {
runtime.primal = basis.recomputex(runtime.instance);
}
Expand Down
17 changes: 0 additions & 17 deletions src/qpsolver/result.hpp
View file

This file was deleted.

20 changes: 8 additions & 12 deletions src/qpsolver/runtime.hpp
View file
Expand Up @@ -6,16 +6,7 @@
#include "instance.hpp"
#include "settings.hpp"
#include "statistics.hpp"

enum class ProblemStatus {
INDETERMINED,
OPTIMAL,
UNBOUNDED,
INFEASIBLE,
ITERATIONLIMIT,
TIMELIMIT,
ERROR
};
#include "qpsolver/qpconst.hpp"

struct Runtime {
Instance instance;
Expand All @@ -32,15 +23,20 @@ struct Runtime {
Vector rowactivity;
Vector dualvar;
Vector dualcon;
ProblemStatus status = ProblemStatus::INDETERMINED;
QpModelStatus status = QpModelStatus::INDETERMINED;

std::vector<BasisStatus> status_var;
std::vector<BasisStatus> status_con;

Runtime(Instance& inst, HighsTimer& ht)
: instance(inst),
timer(ht),
primal(Vector(instance.num_var)),
rowactivity(Vector(instance.num_con)),
dualvar(instance.num_var),
dualcon(instance.num_con) {}
dualcon(instance.num_con),
status_var(instance.num_var),
status_con(instance.num_con) {}
};

#endif

0 comments on commit 975dfc3

Please sign in to comment.