To refine our solutions obtained using heuristics and check out to prove the optimality of solutions, allow us to formulate the Graph Coloring Problem as an ILP model. Remember it may be unable to handle large instances though. The model presented on this section and other exact algorithms are presented in Chapter 3 of Lewis (2021).
Allow us to define the Sets considered on this approach:
- C: Colours
- N: Nodes (or vertices)
- E: Edges
A primary query already comes up “What number of colours should we consider?”. Within the worst-case scenario, all nodes are connected, so one should consider C of the identical size as N. Nonetheless, this approach could make our solutions even harder by unnecessarily increasing the variety of decision variables and worsening the linear leisure of the model. A superb alternative is to make use of a heuristic method (reminiscent of DSatur) to acquire an upper certain for the variety of colours.
On this problem, we now have two groups of decision variables:
- x_{i, c}: Are binary variables indicating that node i is coloured in c
- y_{c}: Are binary variables indicating that color c is used
We must also create constraints to be certain that:
- Every node have to be coloured
- If any node from one edge has a color, be certain that the colour is getting used
- At most 1 node of every edge might be coloured with a given color
- Break symmetry (this will not be required, but might help)
Finally, our objective should minimize the overall variety of colours used, which is the sum of y. Summarizing we now have the next equations.
Without further ado, allow us to import pyomo for the Integer Programming model.
import pyomo.environ as pyo
There are two approaches for modeling an issue in pyomo: Abstract and Concrete models. In the primary approach, the algebraic expressions of the issue are defined before some data values are supplied, whereas, within the second, the model instance is created immediately as its elements are defined. You will discover more about these approaches within the library documentation or within the book by Bynum et al. (2021). Throughout this text, we’ll adopt the Concrete model formulation.
model = pyo.ConcreteModel()
Then, allow us to instantiate our Sets. Parsing the iterables directly from dsatur attributes N and C is valid, so one could use them within the initialize keyword arguments. Alternatively, I’ll pass the unique nodes and edges from our input data and create a variety from the DSatur as our initialization for colours.
model.C = pyo.Set(initialize=range(len(dsatur.C))) # Colours
model.N = pyo.Set(initialize=nodes) # Nodes
model.E = pyo.Set(initialize=edges]) # Edges
Next, we instantiate our decision variables.
model.x = pyo.Var(model.N, model.C, inside=pyo.Binary)
model.y = pyo.Var(model.C, inside=pyo.Binary)
After which our constraints.
# Fill every node with some color
def fill_cstr(model, i):
return sum(model.x[i, :]) == 1# Don't repeat colours on edges and indicator that color is used
def edge_cstr(model, i, j, c):
return model.x[i, c] + model.x[j, c] <= model.y[c]
# Break symmetry by setting a preference order (not required)
def break_symmetry(model, c):
if model.C.first() == c:
return 0 <= model.y[c]
else:
c_prev = model.C.prev(c)
return model.y[c] <= model.y[c_prev]
model.fill_cstr = pyo.Constraint(model.N, rule=fill_cstr)
model.edge_cstr = pyo.Constraint(model.E, model.C, rule=edge_cstr)
model.break_symmetry = pyo.Constraint(model.C, rule=break_symmetry)
You possibly can try to include other symmetry-breaking constraints and see what works higher along with your available solver. In some experiments I performed, including a preference order using total nodes assigned to a given color has been worse than ignoring it. Possibly because of the solver’s native symmetry-breaking strategies.
# Break symmetry by setting a preference order with total assignments
def break_symmetry_agg(model, c):
if model.C.first() == c:
return 0 <= sum(model.x[:, c])
else:
c_prev = model.C.prev(c)
return sum(model.x[:, c]) <= sum(model.x[:, c_prev])
Finally, the target…
# Total variety of colours used
def obj(model):
return sum(model.y[:])model.obj = pyo.Objective(rule=obj)
And our model is able to be solved! To try this I’m using the HiGHS persistent solver, which is offered in pyomo in case highspy can also be installed in your Python environment.
solver = pyo.SolverFactory("appsi_highs")
solver.highs_options["time_limit"] = 120
res = solver.solve(model)
print(res)
For big instances, our solver may need a tough time attempting to improve heuristic solutions. Nonetheless, for a 32-node instance available within the code repository, the solver was capable of reduce the variety of used colours from 9 to eight. Value mentioning that it took 24 seconds to finish execution whereas the DSatur algorithm for a similar instance took only 6 milliseconds (using pure Python, which is an interpreted language).
