From e6f517fca1fe358a48860d7514caae6058b625f6 Mon Sep 17 00:00:00 2001 From: chriseth Date: Wed, 23 Mar 2022 14:54:31 +0100 Subject: [PATCH] delta --- libsolutil/BooleanLP.cpp | 28 +++-- libsolutil/LP.cpp | 260 ++++++++++++++++++++++++++------------- libsolutil/LP.h | 41 ++++-- test/libsolutil/LP.cpp | 6 +- 4 files changed, 229 insertions(+), 106 deletions(-) diff --git a/libsolutil/BooleanLP.cpp b/libsolutil/BooleanLP.cpp index 4577d6582..7b0c97383 100644 --- a/libsolutil/BooleanLP.cpp +++ b/libsolutil/BooleanLP.cpp @@ -233,35 +233,45 @@ pair> BooleanLPSolver::check(vector cons else resizeAndSet(lpState.variableNames, index, name); + // TODO We start afresh here. If we want this to reuse the existing results + // from previous invocations of the boolean solver, we still have to use + // a cache. + // The current optimization is only for CDCL. + m_lpSolver.setState(lpState); + //cout << "Boolean variables:" << joinHumanReadable(booleanVariables) << endl; //cout << "Running LP solver on fixed constraints." << endl; - if (m_lpSolver.check(lpState).first == LPResult::Infeasible) + if (m_lpSolver.check().first == LPResult::Infeasible) { cout << "----->>>>> unsatisfiable" << endl; return {CheckResult::UNSATISFIABLE, {}}; } + set previousConditionalConstraints; auto theorySolver = [&](map const& _booleanAssignment) -> optional { SolvingState lpStateToCheck = lpState; + set conditionalConstraints; for (auto&& [constraintIndex, value]: _booleanAssignment) { if (!value || !state().conditionalConstraints.count(constraintIndex)) continue; + conditionalConstraints.emplace(constraintIndex); + } + set constraintsToRemove = previousConditionalConstraints - conditionalConstraints; + vector constraintsToAdd; + for (size_t constraintIndex: conditionalConstraints - previousConditionalConstraints) + { // "reason" is already stored for those constraints. Constraint const& constraint = state().conditionalConstraints.at(constraintIndex); - solAssert( - constraint.reasons.size() == 1 && - *constraint.reasons.begin() == constraintIndex - ); - lpStateToCheck.constraints.emplace_back(constraint); + solAssert(constraint.reasons.size() == 1 && *constraint.reasons.begin() == constraintIndex); + constraintsToAdd.emplace_back(constraint); } - auto&& [result, modelOrReason] = m_lpSolver.check(move(lpStateToCheck)); + auto&& [result, modelOrReason] = m_lpSolver.check(constraintsToRemove, move(constraintsToAdd)); + previousConditionalConstraints = move(conditionalConstraints); // We can only really use the result "infeasible". Everything else should be "sat". if (result == LPResult::Infeasible) { - // TODO this could be the empty clause if the LP is already infeasible - // with only the fixed constraints - run it beforehand! // TODO is it ok to ignore the non-constraint boolean variables here? Clause conflictClause; for (size_t constraintIndex: get(modelOrReason)) diff --git a/libsolutil/LP.cpp b/libsolutil/LP.cpp index 4ed84c7bb..95f06346b 100644 --- a/libsolutil/LP.cpp +++ b/libsolutil/LP.cpp @@ -22,9 +22,12 @@ #include #include #include +#include + #include #include +#include #include #include #include @@ -502,7 +505,7 @@ bool Constraint::operator<(Constraint const& _other) const if (rational diff = data.get(i) - _other.data.get(i)) return diff < 0; - return false; + return reasons < _other.reasons; } bool Constraint::operator==(Constraint const& _other) const @@ -513,7 +516,8 @@ bool Constraint::operator==(Constraint const& _other) const for (size_t i = 0; i < max(data.size(), _other.data.size()); ++i) if (data.get(i) != _other.data.get(i)) return false; - return true; + + return reasons == _other.reasons; } bool SolvingState::Compare::operator()(SolvingState const& _a, SolvingState const& _b) const @@ -529,6 +533,16 @@ bool SolvingState::Compare::operator()(SolvingState const& _a, SolvingState cons return _a.variableNames < _b.variableNames; } +set SolvingState::reasons() const +{ + set ret; + for (Bounds const& b: bounds) + ret += b.lowerReasons + b.upperReasons; + for (Constraint const& c: constraints) + ret += c.reasons; + return ret; +} + string SolvingState::toString() const { size_t const reasonLength = 10; @@ -719,7 +733,16 @@ void SolvingStateSimplifier::removeEmptyColumns() } } -SolvingState ProblemSplitter::next() +ProblemSplitter::ProblemSplitter(const SolvingState& _state): + m_state(_state), + m_column(1), + m_seenColumns(std::vector(m_state.variableNames.size(), false)) +{ + while (m_column < m_state.variableNames.size() && nonZeroEntriesInColumn(_state, m_column).empty()) + m_column++; +} + +pair, SolvingState> ProblemSplitter::next() { vector includedColumns; vector includedRows; @@ -728,19 +751,12 @@ SolvingState ProblemSplitter::next() // Update state. m_seenColumns |= includedColumns; ++m_column; - while (m_column < m_state.variableNames.size() && m_seenColumns[m_column]) + while (m_column < m_state.variableNames.size() && ( + m_seenColumns[m_column] || + nonZeroEntriesInColumn(m_state, m_column).empty() + )) ++m_column; - if (includedRows.empty()) - { - // This should not happen if the SolvingStateSimplifier has been used beforehand. - // We just check that we did not miss any bounds. - for (auto&& [i, included]: includedColumns | ranges::views::enumerate | ranges::views::tail) - if (included) - solAssert(!m_state.bounds[i].lower && !!m_state.bounds[i].upper); - return next(); - } - SolvingState splitOff; splitOff.variableNames.emplace_back(); @@ -766,100 +782,174 @@ SolvingState ProblemSplitter::next() splitOff.constraints.push_back(move(splitRow)); } - return splitOff; + return {includedColumns, splitOff}; } -LPSolver::LPSolver(bool _supportModels): - m_supportModels(_supportModels), - m_cache(SolvingState::Compare{_supportModels}) +LPSolver::LPSolver(bool) { } -pair> LPSolver::check(SolvingState _state) +void LPSolver::setState(SolvingState _state) { + cout << "Set state:\n" << _state.toString() << endl; + m_subProblems.clear(); + m_subProblemsPerVariable = {}; + m_subProblemsPerConstraint = {}; + normalizeRowLengths(_state); - //cout << "Running LP on:\n" << _state.toString() << endl; - auto&& [simplificationResult, modelOrReasonSet] = SolvingStateSimplifier{_state}.simplify(); - switch (simplificationResult) - { - case LPResult::Infeasible: - //cout << "-> LP infeasible." << endl; - return {LPResult::Infeasible, modelOrReasonSet}; - case LPResult::Feasible: - case LPResult::Unbounded: - solAssert(false); - case LPResult::Unknown: - break; - } + // TODO we should simplify, otherwise we get big problems with constanst that are used everywhere. - Model model = get(modelOrReasonSet); + // TODO we could simplify, but then we need the option to answer 'infeasible' here. + // TODO assert that none of the constraints here have a reason set. - bool canOnlyBeUnknown = false; + cout << "Splitting..." << endl; ProblemSplitter splitter(move(_state)); while (splitter) { - SolvingState split = splitter.next(); - solAssert(!split.constraints.empty(), ""); - solAssert(split.variableNames.size() >= 2, ""); + auto&& [variables, subState] = splitter.next(); + m_subProblems.emplace_back(make_unique(SubProblem{move(subState)})); + for (auto&& [i, included]: variables | ranges::views::enumerate) + if (included) + m_subProblemsPerVariable[i] = m_subProblems.size() - 1; + // We do not need t ofill m_subProblemsPerConstraint because we do not assume + // these constraints to have reasons, so they cannot be removed. + // TODO we cauld assert that + cout << "-------------------- Split out" << endl; + cout << m_subProblems.back()->state.toString() << endl; + } + updateSubProblems(); +} - LPResult lpResult; - vector solution; +pair> LPSolver::check( + std::set const& _constraintsToRemove, + std::vector constraintsToAdd +) +{ + set subproblemsToUpdate; + for (size_t constraintIndex: _constraintsToRemove) + { + cout << "Removing constarint " << constraintIndex << endl; + SubProblem& problem = *m_subProblems[m_subProblemsPerConstraint.at(constraintIndex)]; + problem.dirty = true; + cxx20::erase_if( + problem.state.constraints, + [constraintIndex](Constraint const& _constraint) { + if (_constraint.reasons.count(constraintIndex)) + { + solAssert(_constraint.reasons.size() == 1); + return true; + } + else + return false; + } + ); + } - if (auto conflict = boundsToConstraints(split)) + for (Constraint& constraint: constraintsToAdd) + { + cout << "Adding constraint " << endl; + solAssert(constraint.reasons.size() == 1); + set touchedProblems; + for (auto const& [index, entry]: constraint.data.enumerateTail()) + if (entry && m_subProblemsPerVariable.count(index)) + touchedProblems.emplace(m_subProblemsPerVariable[index]); + if (touchedProblems.empty()) { - //cout << "-> LP infeasible." << endl; - return {LPResult::Infeasible, move(*conflict)}; + m_subProblems.emplace_back(make_unique()); + touchedProblems.emplace(m_subProblems.size() - 1); } + for (size_t problemToErase: touchedProblems | ranges::views::tail | ranges::views::reverse) + combineSubProblems(*touchedProblems.begin(), problemToErase); + addConstraintToSubProblem(*touchedProblems.begin(), move(constraint)); + } + // TODO here, we can try to split again. + // If we split here, then we maybe don't need to split in setState. - auto it = m_cache.find(split); - if (it != m_cache.end()) - tie(lpResult, solution) = it->second; + updateSubProblems(); + + for (unique_ptr const& problem: m_subProblems) + if (problem && problem->result == LPResult::Infeasible) + return {LPResult::Infeasible, problem->state.reasons()}; + + return {LPResult::Unknown, Model{}}; +} + +void LPSolver::combineSubProblems(size_t _combineInto, size_t _combineFrom) +{ + m_subProblems[_combineInto]->dirty = true; + for (auto& item: m_subProblemsPerVariable) + if (item.second == _combineFrom) + item.second = _combineInto; + for (auto& item: m_subProblemsPerConstraint) + if (item.second == _combineFrom) + item.second = _combineInto; + + SolvingState& from = m_subProblems[_combineFrom]->state; + SolvingState& into = m_subProblems[_combineInto]->state; + cout << "Combining " << endl; + cout << from.toString() << "\n------------ with ----------\n" << into.toString() << "\n"; + size_t fromVars = from.variableNames.size(); + size_t intoVars = into.variableNames.size(); + // TODO does it leave out index zero? + into.bounds += from.bounds; + into.variableNames += from.variableNames; + for (Constraint& constraint: into.constraints) + constraint.data.resize(fromVars + intoVars); // -1? + for (Constraint& constraint: from.constraints) + { + Constraint shifted; + shifted.equality = constraint.equality; + shifted.reasons = move(constraint.reasons); + // TODO this can be optimized + for (auto const& [index, entry]: constraint.data.enumerateTail()) + { + shifted.data[index] = {}; + shifted.data[intoVars + index] = entry; + } + into.constraints.emplace_back(move(shifted)); + } + cout << "Result: " << endl; + cout << into.toString() << "\n"; +} + +void LPSolver::addConstraintToSubProblem(size_t _subProblem, Constraint _constraint) +{ + for (auto const& [index, entry]: _constraint.data.enumerateTail()) + if (entry) + { + solAssert(!m_subProblemsPerVariable.count(index) || m_subProblemsPerVariable[index] == _subProblem); + m_subProblemsPerVariable[index] = _subProblem; + } + if (!_constraint.reasons.empty()) + { + solAssert(_constraint.reasons.size() == 1); + m_subProblemsPerConstraint[*_constraint.reasons.begin()] = _subProblem; + } + m_subProblems[_subProblem]->dirty = true; + m_subProblems[_subProblem]->state.constraints.emplace_back(move(_constraint)); +} + +void LPSolver::updateSubProblems() +{ + for (unique_ptr& problem: m_subProblems) + { + if (!problem || !problem->dirty) continue; + + if (auto conflict = boundsToConstraints(problem->state)) + { + problem->result = LPResult::Infeasible; + problem->model = {}; + } else { LinearExpression objectives; objectives.resize(1); - objectives.resize(split.constraints.front().data.size(), rational(bigint(1))); - tie(lpResult, solution) = simplex(split.constraints, move(objectives)); - - // If we do not support models, do not store it in the cache because - // the variable associations will be wrong. - // Otherwise, it is fine to use the model. - m_cache.emplace(split, make_pair(lpResult, m_supportModels ? solution : vector{})); + objectives.resize(problem->state.constraints.front().data.size(), rational(bigint(1))); + tie(problem->result, problem->model) = simplex(problem->state.constraints, move(objectives)); } - - switch (lpResult) - { - case LPResult::Feasible: - case LPResult::Unbounded: + problem->dirty = false; + if (problem->result == LPResult::Infeasible) break; - case LPResult::Infeasible: - { - solAssert(split.bounds.empty()); - set reasons; - for (auto const& constraint: split.constraints) - reasons += constraint.reasons; - //cout << "-> LP infeasible." << endl; - return {LPResult::Infeasible, move(reasons)}; - } - case LPResult::Unknown: - // We do not stop here, because another independent query can still be infeasible. - canOnlyBeUnknown = true; - break; - } - for (auto&& [index, value]: solution | ranges::views::enumerate) - if (index + 1 < split.variableNames.size()) - model[split.variableNames.at(index + 1)] = value; } - - if (canOnlyBeUnknown) - { - //cout << "-> LP unknown." << endl; - return {LPResult::Unknown, Model{}}; - } - - //cout << "-> LP feasible." << endl; - return {LPResult::Feasible, move(model)}; } - - diff --git a/libsolutil/LP.h b/libsolutil/LP.h index bbd0d2cac..7a12f36a9 100644 --- a/libsolutil/LP.h +++ b/libsolutil/LP.h @@ -76,6 +76,8 @@ struct SolvingState // For each bound and constraint, store an index of the literal // that implies it. + std::set reasons() const; + struct Compare { explicit Compare(bool _considerVariableNames = true): considerVariableNames(_considerVariableNames) {} @@ -146,17 +148,13 @@ private: class ProblemSplitter { public: - explicit ProblemSplitter(SolvingState const& _state): - m_state(_state), - m_column(1), - m_seenColumns(std::vector(m_state.variableNames.size(), false)) - {} + explicit ProblemSplitter(SolvingState const& _state); /// @returns true if there are still sub-problems to split out. operator bool() const { return m_column < m_state.variableNames.size(); } /// @returns the next sub-problem. - SolvingState next(); + std::pair, SolvingState> next(); private: SolvingState const& m_state; @@ -182,13 +180,36 @@ class LPSolver public: explicit LPSolver(bool _supportModels = true); - std::pair> check(SolvingState _state); + void setState(SolvingState _state); + /// Modifies the state by removing constraints (identified by their "reason"), + /// adding constraints and then checks for feasibility. + std::pair> check( + std::set const& _constraintsToRemove = {}, + std::vector constraintsToAdd = {} + ); private: - using CacheValue = std::pair>>; + void combineSubProblems(size_t _combineInto, size_t _combineFrom); + void addConstraintToSubProblem(size_t _subProblem, Constraint _constraint); + void updateSubProblems(); - bool m_supportModels = true; - std::map m_cache; + struct SubProblem + { + //std::set variables; + SolvingState state = {}; + bool dirty = true; + LPResult result = LPResult::Unknown; + std::vector> model = {}; + }; + + ReasonSet reasonSetForSubProblem(SubProblem const& _subProblem); + + // TODO we could also use optional + std::vector> m_subProblems; + std::map m_subProblemsPerVariable; + std::map m_subProblemsPerConstraint; + /// TODO also store the first infeasible subproblem? + /// TODO still retain the cache? }; } diff --git a/test/libsolutil/LP.cpp b/test/libsolutil/LP.cpp index b766e9bc7..e9e05de1a 100644 --- a/test/libsolutil/LP.cpp +++ b/test/libsolutil/LP.cpp @@ -92,7 +92,8 @@ public: void feasible(vector> const& _solution) { - auto [result, modelOrReasonSet] = m_solver.check(m_solvingState); + m_solver.setState(m_solvingState); + auto [result, modelOrReasonSet] = m_solver.check(); BOOST_REQUIRE(result == LPResult::Feasible); Model model = get(modelOrReasonSet); for (auto const& [var, value]: _solution) @@ -104,7 +105,8 @@ public: void infeasible(set _reason = {}) { - auto [result, modelOrReason] = m_solver.check(m_solvingState); + m_solver.setState(m_solvingState); + auto [result, modelOrReason] = m_solver.check(); BOOST_REQUIRE(result == LPResult::Infeasible); ReasonSet suppliedReason = get(modelOrReason); BOOST_CHECK_MESSAGE(suppliedReason == _reason, "Reasons are different");