Finding a set of edges to add to a set of graphs to satisfy connectivity constraints

I have a set of N undirected, disconnected graphs that share the same M vertices, but have different edges. I also have a set of constraints for each graph of the form "V₀ is connected to V₁" or "V₃ is not connected to V₅" etc.

I want to find a set of edges such that adding those edges to every graph causes every graph to satisfy its constraints.

For a simple example, consider, for two given graphs with vertices V₀, V₁, V₂, V₃, V₄:

a graph G with one edge (V₀, V₁) and the constraint "V₀ is connected to V₂"
a graph H with edge (V₁, V₃) and the constraints "V₂ is connected to V₄", "V₀ is not connected to V₂"

With these given parameters, the problem is satisfiable by adding the edges {(V₁, V₂), (V₃, V₄)} to both of the graphs.

After failing to solve the problem with a script, I reached for z3 to help, but I've run into trouble trying to express connectivity. My current solution consists of an unquantified formula with only boolean terms:

c_i,j,k for each pair of vertices i, j representing whether they are connected in graph k with added edges from the solution
g_i,j,k for each pair of vertices i, j representing whether the edge between them is given in graph k
s_i,j for each pair of vertices i, j representing whether adding an edge between them is part of the solution

and assertions:

g_i,j,k if the edge (i, j) is given in graph k
¬g_i,j,k if the edge (i, j) is not given in graph k
c_i,j,k if the graph k requires i and j to be connected
¬s_i,j ∧ ¬c_i,j,k if the graph k requires i and j to not be connected
s_i,j ⇒ c_i,j,k for all graphs k
c_i,j,k = s_i,j ∨ g_i,j,k ∨ ((s_i,z ∨ g_i,z,k) ∧ c_z,j,k) for all vertices i, j, z and graphs k

The last clause, however, does not seem to work as intended for expressing connectivity and z3 is returning sat but the clearly necessary s_i,j are false in the resulting model.

Here's the shortest I could manage to make a minimal example of the problem:

from itertools import combinations

from z3 import *

vertices = [0, 1, 2, 3, 4]
edges = [tuple(sorted(edge)) for edge in list(combinations(vertices, 2))]

graphs = {
    "G": [(0, 1)],
    "H": [(1, 3)],
}

facts = Solver()

connected_in_graph = {}
for graph in graphs:
    for edge in edges:
        connected_in_graph[(graph, edge)] = Bool(f"{edge}_conn_in_{graph}")

solution_edges = {}
graph_given_edges = {}
for edge in edges:
    edge_in_solution = Bool(f"{edge}_in_solution")
    solution_edges[edge] = edge_in_solution
    for graph in graphs:
        given = Bool(f"{graph}_{edge}_given")
        graph_given_edges[(graph, edge)] = given
        if edge in graphs[graph]:
            facts.append(given)
        else:
            facts.append(Not(given))

facts.append(
    connected_in_graph[("G", (0, 2))],
    connected_in_graph[("H", (2, 4))],
    Not(connected_in_graph[("H", (0, 2))]),
)

for edge in edges:
    for graph in graphs:
        ors = [
            solution_edges[edge],
            graph_given_edges[(graph, edge)],
        ]
        for vertex in vertices:
            if vertex in edge:
                continue
            l_edge = tuple(sorted((edge[0], vertex)))
            r_edge = tuple(sorted((edge[1], vertex)))
            ors.append(
                And(
                    Or(
                        solution_edges[l_edge],
                        graph_given_edges[(graph, l_edge)],
                    ),
                    connected_in_graph[(graph, r_edge)],
                )
            )
        facts.append(connected_in_graph[(graph, edge)] == Or(*ors))

print(facts.check())
model = facts.model()
for edge in edges:
    c = solution_edges[edge]
    if model[c]:
        print(c)

I suppose what I'd actually need to express, rather than that last constraint, is relations:

c2(i, j, k, through) = s_i,j ∨ g_i,j,k ∨ (∃z. z ∉ (ignored ∪ {i, j}) ∧ (s_i,z ∨ g_i,z,k) ∧ c(z, j, k, ignored ∪ {i}))

c(i, j, k) = c2(i, j, k, {})

However, reducing that to just unquantified boolean terms would obviously take somewhere on the order of M! time and space.

Is there a better way of expressing the existence of a path in a graph in SAT/SMT, or perhaps a better way of solving this problem altogether?

Alias's suggestion to use transitive closure seems to be the solution to this, but I seem to have trouble using it properly. My revised code:

from itertools import combinations

from z3 import *

vertices = [0, 1, 2, 3, 4]
edges = [tuple(sorted(edge)) for edge in list(combinations(vertices, 2))]

graphs = {
    "G": [(0, 1)],
    "H": [(1, 3)],
}

facts = Solver()

vs = {}
VertexSort = DeclareSort("VertexSort")
for vertex in vertices:
    vs[vertex] = Const(f"V{vertex}", VertexSort)
facts.add(Distinct(*vs.values()))

given = {}
directly_connected = {}
tc_connected = {}
for graph in graphs:
    given[graph] = Function(
        f"{graph}_given", VertexSort, VertexSort, BoolSort()
    )
    directly_connected[graph] = Function(
        f"directly_connected_{graph}", VertexSort, VertexSort, BoolSort()
    )
    tc_connected[graph] = TransitiveClosure(directly_connected[graph])

in_solution = Function("in_solution", VertexSort, VertexSort, BoolSort())


for edge in edges:
    # commutativity
    facts.add(
        in_solution(vs[edge[0]], vs[edge[1]])
        == in_solution(vs[edge[1]], vs[edge[0]])
    )
    for graph in graphs:
        # commutativity
        facts.add(
            given[graph](vs[edge[0]], vs[edge[1]])
            == given[graph](vs[edge[1]], vs[edge[0]])
        )
        facts.add(
            directly_connected[graph](vs[edge[0]], vs[edge[1]])
            == directly_connected[graph](vs[edge[1]], vs[edge[0]])
        )
        # definition of direct connection
        facts.add(
            directly_connected[graph](vs[edge[0]], vs[edge[1]])
            == Or(
                in_solution(vs[edge[0]], vs[edge[1]]),
                given[graph](vs[edge[0]], vs[edge[1]]),
            ),
        )
        if edge in graphs[graph]:
            facts.add(given[graph](vs[edge[0]], vs[edge[1]]))
        else:
            facts.add(Not(given[graph](vs[edge[0]], vs[edge[1]])))


facts.append(
    tc_connected["G"](vs[0], vs[2]),
    tc_connected["H"](vs[2], vs[4]),
    Not(tc_connected["H"](vs[0], vs[2])),
)

print(facts.check())
model = facts.model()
print(model)
print(f"solution: {model[in_solution]}")

Prints sat but finds the definition in_solution = [else -> False] rather than the solution I'm looking for. What am I doing wrong?

Solution

As suggested by @alias in this comment, solving the problem was made possible by using transitive closures.

By:

defining relations given_G(v1, v2) and directly_connected_G(v1, v2) for each graph G and a relation in_solution(v1, v2) where v1, v2 are of an enumeration type with constructors for each vertex
asserting directly_connected_G(v1, v2) = in_solution(v1, v2) ∨ given_G(v1, v2) for each v1, v2
declaring a transitive closure TC_connected_G(v1, v2) for each graph G
asserting constraints in terms of TC_connected

I got the solver to return a correct solution for all of my test cases.

The revised code:

from itertools import combinations

from z3 import *

vertices = [0, 1, 2, 3, 4]
edges = [tuple(sorted(edge)) for edge in list(combinations(vertices, 2))]

graphs = {
    "G": [(0, 1)],
    "H": [(1, 3)],
}

facts = Solver()


VertexSort = Datatype("VertexSort")
for vertex in vertices:
    VertexSort.declare(str(vertex))
VertexSort = VertexSort.create()

vs = {}
for vertex in vertices:
    vs[vertex] = getattr(VertexSort, str(vertex))


given = {}
directly_connected = {}
tc_connected = {}
for graph in graphs:
    given[graph] = Function(
        f"{graph}_given", VertexSort, VertexSort, BoolSort()
    )
    directly_connected[graph] = Function(
        f"directly_connected_{graph}", VertexSort, VertexSort, BoolSort()
    )
    tc_connected[graph] = TransitiveClosure(directly_connected[graph])

in_solution = Function("in_solution", VertexSort, VertexSort, BoolSort())


for edge in edges:
    # commutativity
    facts.add(
        in_solution(vs[edge[0]], vs[edge[1]])
        == in_solution(vs[edge[1]], vs[edge[0]])
    )
    for graph in graphs:
        # commutativity
        facts.add(
            given[graph](vs[edge[0]], vs[edge[1]])
            == given[graph](vs[edge[1]], vs[edge[0]])
        )
        facts.add(
            directly_connected[graph](vs[edge[0]], vs[edge[1]])
            == directly_connected[graph](vs[edge[1]], vs[edge[0]])
        )
        # definition of direct connection
        facts.add(
            directly_connected[graph](vs[edge[0]], vs[edge[1]])
            == Or(
                in_solution(vs[edge[0]], vs[edge[1]]),
                given[graph](vs[edge[0]], vs[edge[1]]),
            ),
        )
        if edge in graphs[graph]:
            facts.add(given[graph](vs[edge[0]], vs[edge[1]]))
        else:
            facts.add(Not(given[graph](vs[edge[0]], vs[edge[1]])))


facts.append(
    tc_connected["G"](vs[0], vs[2]),
    tc_connected["H"](vs[2], vs[4]),
    Not(tc_connected["H"](vs[0], vs[2])),
)

print(facts.check())
model = facts.model()
print(f"solution: {model[in_solution]}")