Skip to content

Commit 4317d13

Browse files
levnachlevnach
committed
refactor: move gomory functionaly from int_solver to gomory
1 parent e28b644 commit 4317d13

File tree

4 files changed

+176
-180
lines changed

4 files changed

+176
-180
lines changed

‎src/math/lp/gomory.cpp

+168-15
Original file line numberDiff line numberDiff line change
@@ -401,24 +401,177 @@ class create_cut {
401401
m_f(fractional_part(get_value(basic_inf_int_j).x)),
402402
m_one_minus_f(1 - m_f) {}
403403

404-
};
404+
};
405405

406-
lia_move gomory::cut(lar_term & t, mpq & k, explanation* ex, unsigned basic_inf_int_j, const row_strip<mpq>& row) {
407-
create_cut cc(t, k, ex, basic_inf_int_j, row, lia);
408-
return cc.cut();
409-
}
406+
lia_move gomory::cut(lar_term & t, mpq & k, explanation* ex, unsigned basic_inf_int_j, const row_strip<mpq>& row) {
407+
create_cut cc(t, k, ex, basic_inf_int_j, row, lia);
408+
return cc.cut();
409+
}
410410

411-
lia_move gomory::get_cut(lpvar j) {
412-
unsigned r = lia.row_of_basic_column(j);
413-
const row_strip<mpq>& row = lra.get_row(r);
414-
SASSERT(lra.row_is_correct(r));
415-
SASSERT(lia.is_gomory_cut_target(j));
416-
lia.m_upper = false;
417-
lia.m_cut_vars.push_back(j);
418-
return cut(lia.m_t, lia.m_k, lia.m_ex, j, row);
419-
}
411+
bool gomory::is_gomory_cut_target(lpvar k) {
412+
SASSERT(lia.is_base(k));
413+
// All non base variables must be at their bounds and assigned to rationals (that is, infinitesimals are not allowed).
414+
const row_strip<mpq>& row = lra.get_row(lia.row_of_basic_column(k));
415+
unsigned j;
416+
for (const auto & p : row) {
417+
j = p.var();
418+
if ( k != j && (!lia.at_bound(j) || lia.get_value(j).y != 0)) {
419+
TRACE("gomory_cut", tout << "row is not gomory cut target:\n";
420+
lia.display_column(tout, j);
421+
tout << "infinitesimal: " << !(lia.get_value(j).y ==0) << "\n";);
422+
return false;
423+
}
424+
}
425+
return true;
426+
}
427+
428+
429+
lia_move gomory::get_cut(lpvar j) {
430+
unsigned r = lia.row_of_basic_column(j);
431+
const row_strip<mpq>& row = lra.get_row(r);
432+
SASSERT(lra.row_is_correct(r));
433+
SASSERT(is_gomory_cut_target(j));
434+
lia.m_upper = false;
435+
return cut(lia.m_t, lia.m_k, lia.m_ex, j, row);
436+
}
437+
438+
// return the minimal distance from the variable value to an integer
439+
mpq get_gomory_score(const int_solver& lia, lpvar j) {
440+
const mpq& val = lia.get_value(j).x;
441+
auto l = val - floor(val);
442+
if (l <= mpq(1, 2))
443+
return l;
444+
return mpq(1) - l;
445+
}
446+
447+
unsigned_vector gomory::gomory_select_int_infeasible_vars(unsigned num_cuts) {
448+
std::list<lpvar> sorted_vars;
449+
std::unordered_map<lpvar, mpq> score;
450+
for (lpvar j : lra.r_basis()) {
451+
if (!lia.column_is_int_inf(j) || !is_gomory_cut_target(j))
452+
continue;
453+
SASSERT(!lia.is_fixed(j));
454+
sorted_vars.push_back(j);
455+
score[j] = get_gomory_score(lia, j);
456+
}
457+
// prefer the variables with the values close to integers
458+
sorted_vars.sort([&](lpvar j, lpvar k) {
459+
auto diff = score[j] - score[k];
460+
if (diff.is_neg())
461+
return true;
462+
if (diff.is_pos())
463+
return false;
464+
return lra.usage_in_terms(j) > lra.usage_in_terms(k);
465+
});
466+
unsigned_vector ret;
467+
unsigned n = static_cast<unsigned>(sorted_vars.size());
468+
469+
while (num_cuts-- && n > 0) {
470+
unsigned k = lia.random() % n;
471+
472+
double k_ratio = k / (double) n;
473+
k_ratio *= k_ratio*k_ratio; // square k_ratio to make it smaller
474+
k = static_cast<unsigned>(std::floor(k_ratio * n));
475+
// these operations move k to the beginning of the indices range
476+
SASSERT(0 <= k && k < n);
477+
auto it = sorted_vars.begin();
478+
while(k--) it++;
479+
480+
ret.push_back(*it);
481+
sorted_vars.erase(it);
482+
n--;
483+
}
484+
return ret;
485+
}
486+
487+
488+
lia_move gomory::get_gomory_cuts(unsigned num_cuts) {
489+
struct ex { explanation m_ex; lar_term m_term; mpq m_k; bool m_is_upper; };
490+
unsigned_vector columns_for_cuts = gomory_select_int_infeasible_vars(num_cuts);
491+
492+
vector<ex> cuts;
493+
494+
for (unsigned j : columns_for_cuts) {
495+
lia.m_ex->clear();
496+
lia.m_t.clear();
497+
lia.m_k.reset();
498+
auto r = get_cut(j);
499+
if (r != lia_move::cut)
500+
continue;
501+
cuts.push_back({ *lia.m_ex, lia.m_t, lia.m_k, lia.is_upper() });
502+
if (lia.settings().get_cancel_flag())
503+
return lia_move::undef;
504+
}
420505

506+
auto is_small_cut = [&](ex const& cut) {
507+
return all_of(cut.m_term, [&](auto ci) { return ci.coeff().is_small(); });
508+
};
509+
510+
auto add_cut = [&](ex const& cut) {
511+
u_dependency* dep = nullptr;
512+
for (auto c : cut.m_ex)
513+
dep = lra.join_deps(lra.dep_manager().mk_leaf(c.ci()), dep);
514+
lp::lpvar term_index = lra.add_term(cut.m_term.coeffs_as_vector(), UINT_MAX);
515+
term_index = lra.map_term_index_to_column_index(term_index);
516+
lra.update_column_type_and_bound(term_index,
517+
cut.m_is_upper ? lp::lconstraint_kind::LE : lp::lconstraint_kind::GE,
518+
cut.m_k, dep);
519+
};
520+
521+
auto _check_feasible = [&](void) {
522+
lra.find_feasible_solution();
523+
if (!lra.is_feasible() && !lia.settings().get_cancel_flag()) {
524+
lra.get_infeasibility_explanation(*lia.m_ex);
525+
return false;
526+
}
527+
return true;
528+
};
421529

422-
gomory::gomory(int_solver& lia): lia(lia), lra(lia.lra) { }
530+
bool has_small = false, has_large = false;
531+
532+
for (auto const& cut : cuts) {
533+
if (!is_small_cut(cut)) {
534+
has_large = true;
535+
continue;
536+
}
537+
has_small = true;
538+
add_cut(cut);
539+
}
540+
541+
if (has_large) {
542+
lra.push();
543+
544+
for (auto const& cut : cuts)
545+
if (!is_small_cut(cut))
546+
add_cut(cut);
423547

548+
bool feas = _check_feasible();
549+
lra.pop(1);
550+
551+
if (lia.settings().get_cancel_flag())
552+
return lia_move::undef;
553+
554+
if (!feas)
555+
return lia_move::conflict;
556+
}
557+
558+
if (!_check_feasible())
559+
return lia_move::conflict;
560+
561+
562+
lia.m_ex->clear();
563+
lia.m_t.clear();
564+
lia.m_k.reset();
565+
if (!lia.has_inf_int())
566+
return lia_move::sat;
567+
568+
if (has_small || has_large)
569+
return lia_move::continue_with_check;
570+
571+
lra.move_non_basic_columns_to_bounds();
572+
return lia_move::undef;
573+
}
574+
575+
576+
gomory::gomory(int_solver& lia): lia(lia), lra(lia.lra) { }
424577
}

‎src/math/lp/gomory.h

+4-1
Original file line numberDiff line numberDiff line change
@@ -28,8 +28,11 @@ namespace lp {
2828
class int_solver& lia;
2929
class lar_solver& lra;
3030
lia_move cut(lar_term & t, mpq & k, explanation* ex, unsigned basic_inf_int_j, const row_strip<mpq>& row);
31+
unsigned_vector gomory_select_int_infeasible_vars(unsigned num_cuts);
32+
bool is_gomory_cut_target(lpvar j);
33+
lia_move get_cut(lpvar j);
3134
public:
35+
lia_move gomory::get_gomory_cuts(unsigned num_cuts);
3236
gomory(int_solver& lia);
33-
lia_move get_cut(lpvar j);
3437
};
3538
}

‎src/math/lp/int_solver.cpp

+3-160
Original file line numberDiff line numberDiff line change
@@ -198,8 +198,7 @@ namespace lp {
198198
if (r == lia_move::undef) lra.move_non_basic_columns_to_bounds();
199199
if (r == lia_move::undef && should_hnf_cut()) r = hnf_cut();
200200

201-
std::function<lia_move(lpvar)> gomory_fn = [&](lpvar j) { return gomory(*this).get_cut(j); };
202-
if (r == lia_move::undef && should_gomory_cut()) r = local_cut(2, gomory_fn);
201+
if (r == lia_move::undef && should_gomory_cut()) r = gomory(*this).get_gomory_cuts(2);
203202

204203
if (r == lia_move::undef) r = int_branch(*this)();
205204
if (settings().get_cancel_flag()) r = lia_move::undef;
@@ -692,7 +691,6 @@ namespace lp {
692691
}
693692

694693
void int_solver::simplify(std::function<bool(unsigned)>& is_root) {
695-
696694
return;
697695

698696
#if 0
@@ -829,162 +827,7 @@ namespace lp {
829827

830828
#endif
831829
}
832-
// return the minimal distance from the column value to an integer
833-
mpq get_gomory_score(const int_solver& lia, lpvar j) {
834-
const mpq& val = lia.get_value(j).x;
835-
auto l = val - floor(val);
836-
if (l <= mpq(1, 2))
837-
return l;
838-
return mpq(1) - l;
839-
}
840-
841-
unsigned_vector int_solver::gomory_select_int_infeasible_vars(unsigned num_cuts) {
842-
SASSERT(m_cut_vars.size() == 0&& num_cuts >= 0);
843-
844-
std::list<lpvar> sorted_vars;
845-
std::unordered_map<lpvar, mpq> score;
846-
for (lpvar j : lra.r_basis()) {
847-
if (!column_is_int_inf(j) || !is_gomory_cut_target(j))
848-
continue;
849-
SASSERT(!is_fixed(j));
850-
sorted_vars.push_back(j);
851-
score[j] = get_gomory_score(*this, j);
852-
}
853-
// prefer the columns with the values close to integers
854-
sorted_vars.sort([&](lpvar j, lpvar k) {
855-
auto diff = score[j] - score[k];
856-
if (diff.is_neg())
857-
return true;
858-
if (diff.is_pos())
859-
return false;
860-
return lra.usage_in_terms(j) > lra.usage_in_terms(k);
861-
});
862-
unsigned_vector ret;
863-
unsigned n = static_cast<unsigned>(sorted_vars.size());
864-
865-
while (num_cuts-- && n > 0) {
866-
unsigned k = random() % n;
867-
868-
double k_ratio = k / (double) n;
869-
k_ratio *= k_ratio*k_ratio; // square k_ratio to make it smaller
870-
k = static_cast<unsigned>(std::floor(k_ratio * n));
871-
// these operations move k to the beginning of the indices range
872-
SASSERT(0 <= k && k < n);
873-
auto it = sorted_vars.begin();
874-
while(k--) it++;
875-
876-
ret.push_back(*it);
877-
sorted_vars.erase(it);
878-
n--;
879-
}
880-
return ret;
881-
}
882830

883-
lia_move int_solver::local_cut(unsigned num_cuts, std::function<lia_move(lpvar)>& cut_fn) {
884-
885-
struct ex { explanation m_ex; lar_term m_term; mpq m_k; bool m_is_upper; };
886-
unsigned_vector columns_for_cuts = gomory_select_int_infeasible_vars(num_cuts);
887-
888-
vector<ex> cuts;
889-
890-
for (unsigned j : columns_for_cuts) {
891-
m_ex->clear();
892-
m_t.clear();
893-
m_k.reset();
894-
auto r = cut_fn(j);
895-
if (r != lia_move::cut)
896-
continue;
897-
cuts.push_back({ *m_ex, m_t, m_k, is_upper() });
898-
if (settings().get_cancel_flag())
899-
return lia_move::undef;
900-
}
901-
m_cut_vars.reset();
902-
903-
auto is_small_cut = [&](ex const& cut) {
904-
return all_of(cut.m_term, [&](auto ci) { return ci.coeff().is_small(); });
905-
};
906-
907-
auto add_cut = [&](ex const& cut) {
908-
u_dependency* dep = nullptr;
909-
for (auto c : cut.m_ex)
910-
dep = lra.join_deps(lra.dep_manager().mk_leaf(c.ci()), dep);
911-
lp::lpvar term_index = lra.add_term(cut.m_term.coeffs_as_vector(), UINT_MAX);
912-
term_index = lra.map_term_index_to_column_index(term_index);
913-
lra.update_column_type_and_bound(term_index,
914-
cut.m_is_upper ? lp::lconstraint_kind::LE : lp::lconstraint_kind::GE,
915-
cut.m_k, dep);
916-
};
917-
918-
auto _check_feasible = [&](void) {
919-
lra.find_feasible_solution();
920-
if (!lra.is_feasible() && !settings().get_cancel_flag()) {
921-
lra.get_infeasibility_explanation(*m_ex);
922-
return false;
923-
}
924-
return true;
925-
};
926-
927-
bool has_small = false, has_large = false;
928-
929-
for (auto const& cut : cuts) {
930-
if (!is_small_cut(cut)) {
931-
has_large = true;
932-
continue;
933-
}
934-
has_small = true;
935-
add_cut(cut);
936-
}
937-
938-
if (has_large) {
939-
lra.push();
940-
941-
for (auto const& cut : cuts)
942-
if (!is_small_cut(cut))
943-
add_cut(cut);
944-
945-
bool feas = _check_feasible();
946-
lra.pop(1);
947-
948-
if (settings().get_cancel_flag())
949-
return lia_move::undef;
950-
951-
if (!feas)
952-
return lia_move::conflict;
953-
}
954-
955-
if (!_check_feasible())
956-
return lia_move::conflict;
957-
958-
959-
m_ex->clear();
960-
m_t.clear();
961-
m_k.reset();
962-
if (!has_inf_int())
963-
return lia_move::sat;
964-
965-
if (has_small || has_large)
966-
return lia_move::continue_with_check;
967-
968-
lra.move_non_basic_columns_to_bounds();
969-
return lia_move::undef;
970-
}
971-
972-
bool int_solver::is_gomory_cut_target(lpvar k) {
973-
SASSERT(is_base(k));
974-
// All non base variables must be at their bounds and assigned to rationals (that is, infinitesimals are not allowed).
975-
const row_strip<mpq>& row = lra.get_row(row_of_basic_column(k));
976-
unsigned j;
977-
for (const auto & p : row) {
978-
j = p.var();
979-
if ( k != j && (!at_bound(j) || !is_zero(get_value(j).y))) {
980-
TRACE("gomory_cut", tout << "row is not gomory cut target:\n";
981-
display_column(tout, j);
982-
tout << "infinitesimal: " << !is_zero(get_value(j).y) << "\n";);
983-
return false;
984-
}
985-
}
986-
return true;
987-
}
988-
989-
831+
832+
990833
}

0 commit comments

Comments
 (0)