3
0
Fork 0
mirror of https://github.com/Z3Prover/z3 synced 2025-04-22 16:45:31 +00:00

adding review notes to code

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>
This commit is contained in:
Nikolaj Bjorner 2013-10-26 16:24:21 +08:00
commit 9903c722af
7 changed files with 215 additions and 103 deletions

View file

@ -932,21 +932,23 @@ public:
}
// Return true if there is an edge source --> target.
// If there is such edge, return it in parameter e.
bool get_edge(dl_var source, dl_var target, edge & e) {
// If there is such edge, return its edge_id in parameter id.
bool get_edge_id(dl_var source, dl_var target, edge_id & id) {
edge_id_vector & edges = m_out_edges[source];
typename edge_id_vector::iterator it = edges.begin();
typename edge_id_vector::iterator end = edges.end();
bool found = false;
for (; it != end; ++it) {
edge_id e_id = *it;
edge & e0 = m_edges[e_id];
if (e0.is_enabled() && e0.get_target() == target && !found) {
e = e0;
found = true;
id = *it;
edge & e = m_edges[id];
if (e.is_enabled() && e.get_target() == target) {
return true;
}
}
return found;
return false;
}
edges const & get_all_edges() const {
return m_edges;
}
template<typename Functor>

View file

@ -38,31 +38,30 @@ namespace smt {
// Solve minimum cost flow problem using Network Simplex algorithm
template<typename Ext>
class network_flow : private Ext {
struct GExt : public Ext {
typedef literal explanation;
};
typedef dl_var node;
typedef dl_edge<GExt> edge;
typedef dl_graph<GExt> graph;
typedef typename Ext::numeral numeral;
typedef dl_edge<Ext> edge;
typedef dl_graph<Ext> graph;
typedef typename Ext::numeral numeral;
typedef typename Ext::fin_numeral fin_numeral;
graph m_graph;
// Denote supply/demand b_i on node i
vector<numeral> m_balances;
// Duals of flows which are convenient to compute dual solutions
vector<numeral> m_potentials;
// Keep optimal solution of the min cost flow problem
inf_rational m_objective;
inf_int_rational m_objective;
// Costs on edges
vector<fin_numeral> const & m_costs;
// Basic feasible flows
vector<numeral> m_flows;
// Denote the spanning tree of basic edges
vector<edge> m_basics;
// Denote non-basic edges with flow 0 for uncapicitated networks
vector<edge> m_nonbasics;
// An element is true if the corresponding edge points upwards (compared to the root node)
svector<bool> m_upwards;
// Store the parent of a node in the spanning tree
svector<node> m_pred;
@ -71,7 +70,12 @@ namespace smt {
// Store the pointer to the next node in depth first search ordering
svector<node> m_thread;
bool m_is_optimal;
public:
network_flow(graph & g, vector<fin_numeral> const & costs);
// Initialize the network with a feasible spanning tree
void initialize();
@ -81,13 +85,13 @@ namespace smt {
// If all reduced costs are non-negative, the current flow is optimal
// If not optimal, return a violating edge in the corresponding variable
bool is_optimal(edge & violating_edge);
bool is_optimal(edge_id & violating_edge);
// Send as much flow as possible around the cycle, the first basic edge with flow 0 will leave
edge choose_leaving_edge(const edge & entering_edge);
void update_basics(const edge & entering_edge, const edge & leaving_edge);
edge_id choose_leaving_edge(edge_id entering_edge);
void update_spanning_tree(edge_id entering_edge, edge_id leaving_edge);
bool is_unbounded();
// Compute the optimal solution

View file

@ -24,6 +24,12 @@ Notes:
namespace smt {
template<typename Ext>
network_flow<Ext>::network_flow(graph & g, vector<fin_numeral> const& costs) :
m_graph(g),
m_costs(costs) {
}
template<typename Ext>
void network_flow<Ext>::initialize() {
// TODO: construct an initial spanning tree i.e. inializing m_pred, m_depth and m_thread.
@ -36,54 +42,73 @@ namespace smt {
SASSERT(!m_potentials.empty());
SASSERT(!m_thread.empty());
SASSERT(m_thread.size() == m_pred.size());
array<rational, m_potentials.size()> potentials;
std::copy(m_potentials.begin(), m_potentials.end(), potentials);
rational zero(0);
potentials[0] = zero;
node next = m_thread[0];
while (next != 0) {
node current = m_pred[next];
edge e;
if (m_graph.get_edge(current, next, e)) {
potentials[next] = potentials[current] - e.get_weight();
m_potentials.set(0, numeral::zero());
node target = m_thread[0];
while (target != 0) {
node source = m_pred[target];
edge_id e_id;
if (m_graph.get_edge_id(source, target, e_id)) {
m_potentials.set(target, m_potentials[source] - m_graph.get_weight(e_id));
}
if (m_graph.get_edge(next, current, e)) {
potentials[next] = potentials[current] + e.get_weight();
if (m_graph.get_edge_id(target, source, e_id)) {
m_potentials.set(target, m_potentials[source] + m_graph.get_weight(e_id));
}
next = m_thread[next];
target = m_thread[target];
}
std::copy(potentials.begin(), potentials.end(), m_potentials);
}
template<typename Ext>
void network_flow<Ext>::compute_flows() {
vector<numeral> balances(m_balances);
numeral zero;
m_flows.fill(zero);
vector<edge> basics(m_basics);
// TODO: need a way to find a leaf node of a spanning tree
// OPTIMIZE: Need a set data structure for efficiently removing elements
vector<edge_id> basics;
while (!basics.empty()) {
return;
// Find a leaf node of a spanning tree
node target;
for (unsigned int i = 0; i < m_thread.size(); ++i) {
if (m_depth[i] <= m_depth[m_thread[i]]) {
target = i;
break;
}
}
node source = m_pred[target];
edge_id e_id;
if (m_graph.get_edge_id(source, target, e_id)) {
m_flows.set(e_id, -m_graph.get_weight(basics[target]));
basics[source] += basics[target];
basics.erase(e_id);
}
else if (m_graph.get_edge_id(target, source, e_id)) {
m_flows.set(e_id, m_graph.get_weight(basics[target]));
basics[source] += basics[target];
basics.erase(e_id);
}
}
}
template<typename Ext>
bool network_flow<Ext>::is_optimal(edge & violating_edge) {
typename vector<edge>::iterator it = m_nonbasics.begin();
typename vector<edge>::iterator end = m_nonbasics.end();
bool network_flow<Ext>::is_optimal(edge_id & violating_edge) {
// TODO: how to get nonbasics vector?
vector<edge> nonbasics;
typename vector<edge>::iterator it = nonbasics.begin();
typename vector<edge>::iterator end = nonbasics.end();
bool found = false;
for (unsigned int i = 0; i < m_nonbasics.size(); ++i) {
edge & e = m_nonbasics[i];
for (unsigned int i = 0; i < nonbasics.size(); ++i) {
edge & e = nonbasics[i];
if (e.is_enabled()) {
node source = e.get_source();
node target = e.get_target();
numeral cost = e.get_weight() - m_potentials[source] + m_potentials[target];
// Choose the first negative-cost edge to be the violating edge
// TODO: add multiple pivoting strategies
if (cost < 0) {
violating_edge = e;
numeral zero(0);
if (cost < zero) {
edge_id e_id;
m_graph.get_edge_id(source, target, e_id);
violating_edge = e_id;
found = true;
break;
}
@ -93,9 +118,9 @@ namespace smt {
}
template<typename Ext>
dl_edge<typename network_flow<Ext>::GExt> network_flow<Ext>::choose_leaving_edge(const edge & entering_edge) {
node source = entering_edge.get_source();
node target = entering_edge.get_target();
edge_id network_flow<Ext>::choose_leaving_edge(edge_id entering_edge) {
node source = m_graph.get_source(entering_edge);
node target = m_graph.get_target(entering_edge);
while (source != target) {
if (m_depth[source] > m_depth[target])
source = m_pred[source];
@ -106,14 +131,28 @@ namespace smt {
target = m_pred[target];
}
}
edge e;
m_graph.get_edge(source, target, e);
return e;
edge_id e_id;
m_graph.get_edge_id(source, target, e_id);
return e_id;
}
template<typename Ext>
void network_flow<Ext>::update_basics(const edge & entering_edge, const edge & leaving_edge) {
void network_flow<Ext>::update_spanning_tree(edge_id entering_edge, edge_id leaving_edge) {
// Need special handling in case two edges are identical
SASSERT(entering_edge != leaving_edge);
// Update potentials
node target = m_upwards[leaving_edge] ? m_graph.get_source(leaving_edge) : m_graph.get_target(leaving_edge);
numeral src_pot = m_potentials[m_graph.get_source(entering_edge)];
numeral tgt_pot = m_potentials[m_graph.get_target(entering_edge)];
numeral weight = m_graph.get_weight(entering_edge);
numeral change = m_upwards[entering_edge] ? (weight - src_pot + tgt_pot) : (-weight + src_pot - tgt_pot);
m_potentials[target] += change;
node start = m_thread[target];
while (m_depth[start] > m_depth[target]) {
m_potentials[start] += change;
start = m_thread[start];
}
}
template<typename Ext>
@ -124,24 +163,34 @@ namespace smt {
// Get the optimal solution
template<typename Ext>
void network_flow<Ext>::get_optimal_solution(numeral & objective, vector<numeral> & flows) {
SASSERT(m_is_optimal);
flows.reset();
flows.append(m_flows);
// TODO: calculate objective value
numeral cost(0);
for (unsigned int i = 0; i < m_flows.size(); ++i) {
// FIXME: this * operator is not supported
cost += m_costs[i] * m_flows[i];
}
objective = cost;
}
// Minimize cost flows
// Return true if found an optimal solution, and return false if unbounded
template<typename Ext>
bool network_flow<Ext>::min_cost() {
SASSERT(!m_graph.get_all_edges().empty());
initialize();
edge & entering_edge;
edge_id entering_edge;
while (!is_optimal(entering_edge)) {
edge & leaving_edge = choose_leaving_edge();
update_tree(entering_edge, leaving_edge);
if (is_unbounded())
return false;
edge_id leaving_edge = choose_leaving_edge(entering_edge);
update_spanning_tree(entering_edge, leaving_edge);
if (is_unbounded()) {
m_is_optimal = false;
return m_is_optimal;
}
}
return true;
m_is_optimal = true;
return m_is_optimal;
}
}

View file

@ -38,7 +38,7 @@ Revision History:
#include"numeral_factory.h"
#include"smt_clause.h"
#include"theory_opt.h"
#include"network_flow_def.h"
#include"network_flow.h"
// The DL theory can represent term such as n + k, where n is an enode and k is a numeral.
namespace smt {
@ -307,14 +307,13 @@ namespace smt {
virtual bool maximize(theory_var v);
virtual theory_var add_objective(app* term);
virtual inf_eps_rational<inf_rational> get_objective_value(theory_var v);
numeral m_objective_value;
typedef vector <std::pair<theory_var, rational> > objective_term;
vector<objective_term> m_objectives;
vector<rational> m_objective_vars;
vector<rational> m_objective_consts;
void internalize_objective(app * n, objective_term & objective);
network_flow<Ext> m_network_flow;
bool internalize_objective(expr * n, rational const& m, rational& r, objective_term & objective);
private:
@ -368,6 +367,7 @@ namespace smt {
// TODO: It doesn't need to be a rational, but a bignum integer.
static const bool m_int_theory = true;
typedef rational numeral;
typedef rational fin_numeral;
numeral m_epsilon;
idl_ext() : m_epsilon(1) {}
};
@ -376,6 +376,7 @@ namespace smt {
// TODO: It doesn't need to be a rational, but a bignum integer.
static const bool m_int_theory = true;
typedef s_integer numeral;
typedef s_integer fin_numeral;
numeral m_epsilon;
sidl_ext() : m_epsilon(1) {}
};
@ -383,13 +384,15 @@ namespace smt {
struct rdl_ext {
static const bool m_int_theory = false;
typedef inf_int_rational numeral;
numeral m_epsilon;
typedef rational fin_numeral;
numeral m_epsilon;
rdl_ext() : m_epsilon(rational(), true) {}
};
struct srdl_ext {
static const bool m_int_theory = false;
typedef inf_s_integer numeral;
typedef s_integer fin_numeral;
numeral m_epsilon;
srdl_ext() : m_epsilon(s_integer(0),true) {}
};

View file

@ -28,6 +28,7 @@ Revision History:
#include"ast_pp.h"
#include"warning.h"
#include"smt_model_generator.h"
#include"network_flow_def.h"
using namespace smt;
@ -55,8 +56,6 @@ std::ostream& theory_diff_logic<Ext>::atom::display(theory_diff_logic const& th,
template<typename Ext>
void theory_diff_logic<Ext>::nc_functor::reset() {
m_antecedents.reset();
m_objectives.reset();
m_objective_vars.reset();
}
@ -777,6 +776,8 @@ void theory_diff_logic<Ext>::reset_eh() {
m_agility = 0.5;
m_is_lia = true;
m_non_diff_logic_exprs = false;
m_objectives .reset();
m_objective_consts.reset();
theory::reset_eh();
}
@ -1000,14 +1001,51 @@ void theory_diff_logic<Ext>::get_implied_bound_antecedents(edge_id bridge_edge,
template<typename Ext>
bool theory_diff_logic<Ext>::maximize(theory_var v) {
objective_term const& objective = m_objectives[v];
IF_VERBOSE(1,
objective_term const& objective = m_objectives[v];
for (unsigned i = 0; i < objective.size(); ++i) {
verbose_stream() << "coefficient " << objective[i].second << " of theory_var " << objective[i].first << "\n";
}
verbose_stream() << m_objective_vars[v] << "\n";);
verbose_stream() << m_objective_consts[v] << "\n";);
NOT_IMPLEMENTED_YET();
return false;
// Double the number of edges in the new graph
// NSB review: what about disabled edges? They should not be added, right?
// For efficiency we probably want to reuse m_graph and keep extra edges on the side or add them to
// m_graph as well.
dl_graph<GExt> g;
vector<dl_edge<GExt> > const& es = m_graph.get_all_edges();
dl_var offset = m_graph.get_num_edges();
for (unsigned i = 0; i < es.size(); ++i) {
dl_edge<GExt> const & e = es[i];
if (e.is_enabled()) {
g.enable_edge(g.add_edge(e.get_source(), e.get_target(), e.get_weight(), e.get_explanation()));
g.enable_edge(g.add_edge(e.get_target() + offset, e.get_source() + offset, e.get_weight(), e.get_explanation()));
}
}
// Objective coefficients now become costs
vector<fin_numeral> base_costs, aux_costs;
for (unsigned i = 0; i < objective.size(); ++i) {
fin_numeral cost(objective[i].second);
base_costs.push_back(cost);
aux_costs.push_back(-cost);
}
vector<fin_numeral> costs;
costs.append(base_costs);
costs.append(aux_costs);
network_flow<GExt> net_flow(g, costs);
bool is_optimal = net_flow.min_cost();
if (is_optimal) {
numeral objective_value;
vector<numeral> flows;
net_flow.get_optimal_solution(objective_value, flows);
m_objective_value = objective_value.get_rational();
// TODO: return the model of the optimal solution
}
return is_optimal;
}
template<typename Ext>
@ -1015,59 +1053,55 @@ theory_var theory_diff_logic<Ext>::add_objective(app* term) {
objective_term objective;
theory_var result = m_objectives.size();
rational q(1), r(0);
if (!internalize_objective(term, q, r, objective)) {
result = null_theory_var;
if (internalize_objective(term, q, r, objective)) {
m_objectives.push_back(objective);
m_objective_consts.push_back(r);
}
else {
m_objectives.push_back(objective);
m_objective_vars.push_back(r);
result = null_theory_var;
}
return result;
}
template<typename Ext>
inf_eps_rational<inf_rational> theory_diff_logic<Ext>::get_objective_value(theory_var v) {
NOT_IMPLEMENTED_YET();
inf_rational objective;
inf_eps_rational<inf_rational> val(objective);
return val;
}
template<typename Ext>
bool theory_diff_logic<Ext>::internalize_objective(app * n, rational& q, objective_term & objective) {
bool theory_diff_logic<Ext>::internalize_objective(expr * n, rational const& m, rational& q, objective_term & objective) {
// Compile term into objective_term format
rational r;
expr* x, *y;
if (m_util.is_numeral(n, r)) {
q += r;
}
else if (m_util.is_add(n)) {
for (unsigned i = 0; i < n->get_num_args(); ++i) {
if (!internalize_objective(to_app(n->get_arg(i)), objective)) {
for (unsigned i = 0; i < to_app(n)->get_num_args(); ++i) {
if (!internalize_objective(to_app(n)->get_arg(i), m, q, objective)) {
return false;
}
};
}
}
else if (m_util.is_mul(n)) {
SASSERT(n->get_num_args() == 2);
rational r;
app * lhs = to_app(n->get_arg(0));
app * rhs = to_app(n->get_arg(1));
SASSERT(m_util.is_numeral(lhs, r) || m_util.is_numeral(rhs, r));
if (!m_util.is_numeral(lhs, r))
std::swap(lhs, rhs);
m_util.is_numeral(lhs, r);
theory_var v = mk_var(rhs);
objective.push_back(std::make_pair(v, r));
else if (m_util.is_mul(n, x, y) && m_util.is_numeral(x, r)) {
return internalize_objective(y, m*r, q, objective);
}
else if (n->get_family_id() == m_util.get_family_id()) {
else if (m_util.is_mul(n, y, x) && m_util.is_numeral(x, r)) {
return internalize_objective(y, m*r, q, objective);
}
else if (!is_app(n)) {
return false;
}
else if (to_app(n)->get_family_id() == m_util.get_family_id()) {
return false;
}
else {
theory_var v = mk_var(n);
rational r(1);
objective.push_back(std::make_pair(v, r));
theory_var v = mk_var(to_app(n));
objective.push_back(std::make_pair(v, m));
}
return true;
}

View file

@ -155,6 +155,17 @@ class inf_int_rational {
return *this;
}
inf_int_rational & operator*=(const rational & r) {
if (!r.is_int32()) {
throw default_exception("multiplication with large rational is not possible");
}
m_first *= r;
m_second *= r.get_int32();
return *this;
}
inf_int_rational & operator-=(const inf_int_rational & r) {
m_first -= r.m_first;
m_second -= r.m_second;
@ -344,6 +355,10 @@ inline inf_int_rational operator+(const inf_int_rational & r1, const inf_int_rat
return inf_int_rational(r1) += r2;
}
inline inf_int_rational operator*(const rational & r1, const inf_int_rational & r2) {
return inf_int_rational(r2) *= r1;
}
inline inf_int_rational operator-(const inf_int_rational & r1, const inf_int_rational & r2) {
return inf_int_rational(r1) -= r2;
}

View file

@ -111,6 +111,11 @@ public:
return INT_MIN <= v && v <= INT_MAX;
}
int get_int32() const {
SASSERT(is_int32());
return (int)get_int64();
}
double get_double() const { return m().get_double(m_val); }
rational const & get_rational() const { return *this; }