delta

libsolutil/BooleanLP.cpp

@@ -233,35 +233,45 @@ pair<CheckResult, vector<string>> BooleanLPSolver::check(vector<Expression> 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<size_t> previousConditionalConstraints;
 	auto theorySolver = [&](map<size_t, bool> const& _booleanAssignment) -> optional<Clause>
 	{
 		SolvingState lpStateToCheck = lpState;
+		set<size_t> conditionalConstraints;
 		for (auto&& [constraintIndex, value]: _booleanAssignment)
 		{
 			if (!value || !state().conditionalConstraints.count(constraintIndex))
 				continue;
+			conditionalConstraints.emplace(constraintIndex);
+		}
+		set<size_t> constraintsToRemove = previousConditionalConstraints - conditionalConstraints;
+		vector<Constraint> constraintsToAdd;
+		for (size_t constraintIndex: conditionalConstraints - previousConditionalConstraints)
+		{
 			// "reason" is already stored for those constraints.
 			Constraint const& constraint = state().conditionalConstraints.at(constraintIndex);
-			solAssert(
+			solAssert(constraint.reasons.size() == 1 && *constraint.reasons.begin() == constraintIndex);
-				constraint.reasons.size() == 1 &&
+			constraintsToAdd.emplace_back(constraint);
-				*constraint.reasons.begin() == constraintIndex
-			);
-			lpStateToCheck.constraints.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<ReasonSet>(modelOrReason))

libsolutil/LP.cpp

@@ -22,9 +22,12 @@
 #include <libsolutil/CommonIO.h>
 #include <libsolutil/StringUtils.h>
 #include <libsolutil/LinearExpression.h>
+#include <libsolutil/cxx20.h>
+
 #include <liblangutil/Exceptions.h>
 
 #include <range/v3/view/enumerate.hpp>
+#include <range/v3/view/reverse.hpp>
 #include <range/v3/view/transform.hpp>
 #include <range/v3/view/filter.hpp>
 #include <range/v3/view/tail.hpp>
@@ -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<size_t> SolvingState::reasons() const
+{
+	set<size_t> 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<bool>(m_state.variableNames.size(), false))
+{
+	while (m_column < m_state.variableNames.size() && nonZeroEntriesInColumn(_state, m_column).empty())
+		m_column++;
+}
+
+pair<vector<bool>, SolvingState> ProblemSplitter::next()
 {
 	vector<bool> includedColumns;
 	vector<bool> 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):
+LPSolver::LPSolver(bool)
-	m_supportModels(_supportModels),
-	m_cache(SolvingState::Compare{_supportModels})
 {
 }
 
-pair<LPResult, variant<Model, ReasonSet>> 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();
+	// TODO we should simplify, otherwise we get big problems with constanst that are used everywhere.
-	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;
-	}
 
-	Model model = get<Model>(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();
+		auto&& [variables, subState] = splitter.next();
-		solAssert(!split.constraints.empty(), "");
+		m_subProblems.emplace_back(make_unique<SubProblem>(SubProblem{move(subState)}));
-		solAssert(split.variableNames.size() >= 2, "");
+		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;
+pair<LPResult, variant<Model, ReasonSet>> LPSolver::check(
-		vector<rational> solution;
+	std::set<size_t> const& _constraintsToRemove,
+	std::vector<Constraint> constraintsToAdd
+)
+{
+	set<size_t> 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<size_t> 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;
+			m_subProblems.emplace_back(make_unique<SubProblem>());
-			return {LPResult::Infeasible, move(*conflict)};
+			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);
+	updateSubProblems();
-		if (it != m_cache.end())
+
-			tie(lpResult, solution) = it->second;
+	for (unique_ptr<SubProblem> 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<SubProblem>& 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)));
+			objectives.resize(problem->state.constraints.front().data.size(), rational(bigint(1)));
-			tie(lpResult, solution) = simplex(split.constraints, move(objectives));
+			tie(problem->result, problem->model) = simplex(problem->state.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<rational>{}));
 		}
-
+		problem->dirty = false;
-		switch (lpResult)
+		if (problem->result == LPResult::Infeasible)
-		{
-		case LPResult::Feasible:
-		case LPResult::Unbounded:
 			break;
-		case LPResult::Infeasible:
-		{
-			solAssert(split.bounds.empty());
-			set<size_t> 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)};
 }
-
-

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<size_t> reasons() const;
+
 	struct Compare
 	{
 		explicit Compare(bool _considerVariableNames = true): considerVariableNames(_considerVariableNames) {}
@@ -146,17 +148,13 @@ private:
 class ProblemSplitter
 {
 public:
-	explicit ProblemSplitter(SolvingState const& _state):
+	explicit ProblemSplitter(SolvingState const& _state);
-		m_state(_state),
-		m_column(1),
-		m_seenColumns(std::vector<bool>(m_state.variableNames.size(), false))
-	{}
 
 	/// @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<std::vector<bool>, SolvingState> next();
 
 private:
 	SolvingState const& m_state;
@@ -182,13 +180,36 @@ class LPSolver
 public:
 	explicit LPSolver(bool _supportModels = true);
 
-	std::pair<LPResult, std::variant<Model, ReasonSet>> 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<LPResult, std::variant<Model, ReasonSet>> check(
+		std::set<size_t> const& _constraintsToRemove = {},
+		std::vector<Constraint> constraintsToAdd = {}
+	);
 
 private:
-	using CacheValue = std::pair<LPResult, std::vector<boost::rational<bigint>>>;
+	void combineSubProblems(size_t _combineInto, size_t _combineFrom);
+	void addConstraintToSubProblem(size_t _subProblem, Constraint _constraint);
+	void updateSubProblems();
 
-	bool m_supportModels = true;
+	struct SubProblem
-	std::map<SolvingState, CacheValue, SolvingState::Compare> m_cache;
+	{
+		//std::set<size_t> variables;
+		SolvingState state = {};
+		bool dirty = true;
+		LPResult result = LPResult::Unknown;
+		std::vector<boost::rational<bigint>> model = {};
+	};
+
+	ReasonSet reasonSetForSubProblem(SubProblem const& _subProblem);
+
+	// TODO we could also use optional
+	std::vector<std::unique_ptr<SubProblem>> m_subProblems;
+	std::map<size_t, size_t> m_subProblemsPerVariable;
+	std::map<size_t, size_t> m_subProblemsPerConstraint;
+	/// TODO also store the first infeasible subproblem?
+	/// TODO still retain the cache?
 };
 
 }

test/libsolutil/LP.cpp

@@ -92,7 +92,8 @@ public:
 
 	void feasible(vector<pair<string, rational>> 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<Model>(modelOrReasonSet);
 		for (auto const& [var, value]: _solution)
@@ -104,7 +105,8 @@ public:
 
 	void infeasible(set<size_t> _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<ReasonSet>(modelOrReason);
 		BOOST_CHECK_MESSAGE(suppliedReason == _reason, "Reasons are different");