14 #ifndef LIMBO_ALGORITHMS_COLORING_LP
15 #define LIMBO_ALGORITHMS_COLORING_LP
27 #include <boost/cstdint.hpp>
28 #include <boost/graph/graph_concepts.hpp>
29 #include <boost/property_map/property_map.hpp>
30 #include <boost/dynamic_bitset.hpp>
33 #include "gurobi_c++.h"
35 #ifdef DEBUG_NONINTEGERS
36 extern std::vector<unsigned int> vLP1NonInteger;
37 extern std::vector<unsigned int> vLP1HalfInteger;
38 extern std::vector<unsigned int> vLP2NonInteger;
39 extern std::vector<unsigned int> vLP2HalfInteger;
40 extern std::vector<unsigned int> vLPEndNonInteger;
41 extern std::vector<unsigned int> vLPEndHalfInteger;
42 extern std::vector<unsigned int> vLPNumIter;
66 template <
typename GraphType>
72 typedef typename base_type::graph_type graph_type;
73 typedef typename base_type::graph_vertex_type graph_vertex_type;
74 typedef typename base_type::graph_edge_type graph_edge_type;
75 typedef typename base_type::vertex_iterator_type vertex_iterator_type;
76 typedef typename base_type::edge_iterator_type edge_iterator_type;
77 typedef typename base_type::edge_weight_type edge_weight_type;
91 : vertex_non_integer_num(std::numeric_limits<uint32_t>::
max())
92 , edge_non_integer_num(std::numeric_limits<uint32_t>::
max())
93 , vertex_half_integer_num(std::numeric_limits<uint32_t>::
max())
94 , edge_half_integer_num(std::numeric_limits<uint32_t>::
max())
111 void set(
double c,
char s)
121 if (coeff == 0.0 || rhs.
coeff == 0.0)
123 else if (sense == rhs.
sense)
124 return (coeff > 0 && rhs.
coeff > 0) || (coeff < 0 && rhs.
coeff < 0);
126 return (coeff > 0 && rhs.
coeff < 0) || (coeff < 0 && rhs.
coeff > 0);
161 void set_optimize_model(vector<GRBVar>& vColorBits, vector<GRBVar>& vEdgeBits, GRBLinExpr& obj, GRBModel& optModel);
164 void set_anchor(vector<GRBVar>& vColorBits)
const;
184 void get_odd_cycles(graph_vertex_type
const& v, vector<vector<graph_vertex_type> >& vOddCyle);
204 void non_integer_num(vector<GRBVar>
const& vColorBits, vector<GRBVar>
const& vEdgeBits, NonIntegerInfo& info)
const;
209 void non_integer_num(vector<GRBVar>
const& vVariables, uint32_t& nonIntegerNum, uint32_t& halfIntegerNum)
const;
214 bool is_integer(
double value)
const {
return value == floor(value);}
219 uint32_t
check_precolored_num(vector<LPColoring<GraphType>::graph_vertex_type>
const& vVertex)
const;
228 template <
typename GraphType>
231 , m_vVertexHandledByOddCycle(
boost::num_vertices(g), false)
238 template <
typename GraphType>
245 if (m_vVertexHandledByOddCycle[v])
return;
246 else m_vVertexHandledByOddCycle[v] =
true;
249 uint32_t numVertices = boost::num_vertices(this->m_graph);
250 vector<int32_t> vNodeDistColor (numVertices, -1);
251 vector<bool> vNodeVisited (numVertices,
false);
254 vector<bool> vInCycle (numVertices,
false);
258 vector<graph_vertex_type> vVertexStack;
259 vVertexStack.reserve(numVertices);
260 vNodeVisited[v] =
true;
261 vNodeDistColor[v] = 0;
262 vVertexStack.push_back(v);
263 while (!vVertexStack.empty())
269 graph_vertex_type cv = vVertexStack.back();
271 typename boost::graph_traits<graph_type>::adjacency_iterator ui, uie;
272 for(boost::tie(ui, uie) = adjacent_vertices(cv, this->m_graph); ui != uie; ++ui)
274 graph_vertex_type u = *ui;
276 if(vNodeDistColor[u] == -1)
278 vNodeDistColor[u] = 1 - vNodeDistColor[cv];
279 vNodeVisited[u] =
true;
281 vVertexStack.push_back(u);
289 vector<graph_vertex_type> cycle;
291 for (boost::tie(ui, uie) = adjacent_vertices(cv, this->m_graph); ui != uie; ++ui)
293 graph_vertex_type u = *ui;
294 if (vNodeVisited[u] && (vNodeDistColor[u] == vNodeDistColor[cv]))
300 int32_t cnt = vVertexStack.size();
301 for(int32_t k = cnt - 1; k >= 0; --k)
303 cycle.push_back(vVertexStack[k]);
305 vInCycle[vVertexStack[k]] =
true;
306 if(vVertexStack[k] == u)
break;
310 if (!cycle.empty() && vInCycle[v] && !(this->has_precolored() && check_precolored_num(cycle)+2 > cycle.size()))
312 vOddCyle.push_back(vector<graph_vertex_type>());
313 vOddCyle.back().swap(cycle);
316 for (int32_t k = 0; k != cnt; ++k)
317 m_vVertexHandledByOddCycle[vVertexStack[k]] =
true;
324 vVertexStack.pop_back();
325 vNodeVisited[cv] =
false;
331 template <
typename GraphType>
334 graph_vertex_type u = 0;
335 uint32_t maxDegree = 0;
336 vertex_iterator_type vi, vie;
337 for (boost::tie(vi, vie) = boost::vertices(this->m_graph); vi != vie; ++vi)
339 uint32_t d = boost::degree(*vi, this->m_graph);
349 template <
typename GraphType>
352 uint32_t numVertices = boost::num_vertices(this->m_graph);
354 uint32_t numColorBits = numVertices<<1;
357 vColorBits.reserve(numColorBits);
361 for (uint32_t i = 0; i != numColorBits; ++i)
363 sprintf(buf,
"v%u", i);
365 uint32_t vertexIdx = i>>1;
366 int8_t color = this->m_vColor[vertexIdx];
368 vColorBits.push_back(optModel.addVar(0.0, 1.0, 0.0, GRB_CONTINUOUS, buf));
371 int8_t colorBit = (i&1)? (color&1) : (color>>1);
372 #ifdef DEBUG_LPCOLORING
373 assert(colorBit >= 0 && colorBit <= 1);
375 vColorBits.push_back(optModel.addVar(colorBit, colorBit, colorBit, GRB_CONTINUOUS, buf));
389 optModel.setObjective(obj, GRB_MINIMIZE);
392 edge_iterator_type ei, eie;
393 for(boost::tie(ei, eie) = boost::edges(this->m_graph); ei != eie; ++ei)
395 graph_edge_type e = *ei;
396 graph_vertex_type s = boost::source(e, this->m_graph);
397 graph_vertex_type t = boost::target(e, this->m_graph);
401 if (this->has_precolored() && this->m_vColor[s] >= 0 && this->m_vColor[t] >= 0)
404 uint32_t bitIdxS = s<<1;
405 uint32_t bitIdxT = t<<1;
407 edge_weight_type w = this->edge_weight(e);
408 assert_msg(w > 0,
"no stitch edge allowed, positive edge weight expected: " << w);
410 sprintf(buf,
"R%u", m_constrs_num++);
412 vColorBits[bitIdxS] + vColorBits[bitIdxS+1]
413 + vColorBits[bitIdxT] + vColorBits[bitIdxT+1] >= 1
416 sprintf(buf,
"R%u", m_constrs_num++);
418 1 - vColorBits[bitIdxS] + vColorBits[bitIdxS+1]
419 + 1 - vColorBits[bitIdxT] + vColorBits[bitIdxT+1] >= 1
422 sprintf(buf,
"R%u", m_constrs_num++);
424 vColorBits[bitIdxS] + 1 - vColorBits[bitIdxS+1]
425 + vColorBits[bitIdxT] + 1 - vColorBits[bitIdxT+1] >= 1
428 sprintf(buf,
"R%u", m_constrs_num++);
430 1 - vColorBits[bitIdxS] + 1 - vColorBits[bitIdxS+1]
431 + 1 - vColorBits[bitIdxT] + 1 - vColorBits[bitIdxT+1] >= 1
435 if (this->color_num() == base_type::THREE)
437 for(uint32_t k = 0; k < numColorBits; k += 2)
439 sprintf(buf,
"R%u", m_constrs_num++);
441 vColorBits[k] + vColorBits[k+1] <= 1,
450 template <
typename GraphType>
453 if (this->has_precolored())
456 graph_vertex_type anchorVertex = max_degree_vertex();
457 uint32_t bitIdxAnchor = anchorVertex<<1;
458 vColorBits[bitIdxAnchor].set(GRB_DoubleAttr_UB, 0.0);
459 vColorBits[bitIdxAnchor].set(GRB_DoubleAttr_LB, 0.0);
460 vColorBits[bitIdxAnchor+1].set(GRB_DoubleAttr_UB, 0.0);
461 vColorBits[bitIdxAnchor+1].set(GRB_DoubleAttr_LB, 0.0);
465 template <
typename GraphType>
468 for(uint32_t k = 0, ke = vColorBits.size(); k < ke; k += 2)
470 double value1 = vColorBits[k].get(GRB_DoubleAttr_X);
471 double value2 = vColorBits[k+1].get(GRB_DoubleAttr_X);
475 obj += vColorBits[k+1]-vColorBits[k];
476 else if (value1 < value2)
477 obj += vColorBits[k]-vColorBits[k+1];
483 template <
typename GraphType>
486 edge_iterator_type ei, eie;
487 for(boost::tie(ei, eie) = boost::edges(this->m_graph); ei != eie; ++ei)
489 graph_edge_type e = *ei;
490 graph_vertex_type s = boost::source(e, this->m_graph);
491 graph_vertex_type t = boost::target(e, this->m_graph);
493 for (uint32_t i = 0; i != 2; ++i)
495 uint32_t bitIdxS = (s<<1)+i;
496 uint32_t bitIdxT = (t<<1)+i;
498 double value1 = vColorBits[bitIdxS].get(GRB_DoubleAttr_X);
499 double value2 = vColorBits[bitIdxT].get(GRB_DoubleAttr_X);
502 obj += vColorBits[bitIdxT]-vColorBits[bitIdxS];
503 else if (value1 < value2)
504 obj += vColorBits[bitIdxS]-vColorBits[bitIdxT];
510 template <
typename GraphType>
514 vector<vector<graph_vertex_type> > vOddCyle;
515 for(uint32_t k = 0, ke = vColorBits.size(); k < ke; k += 2)
518 if (vColorBits[k].
get(GRB_DoubleAttr_X) != 0.5 && vColorBits[k+1].get(GRB_DoubleAttr_X) != 0.5)
521 graph_vertex_type v = k>>1;
522 this->get_odd_cycles(v, vOddCyle);
524 for (
typename vector<vector<graph_vertex_type> >::const_iterator it1 = vOddCyle.begin(), it1e = vOddCyle.end(); it1 != it1e; ++it1)
526 vector<graph_vertex_type>
const& oddCycle = *it1;
527 int32_t cycleLength = oddCycle.size();
528 GRBLinExpr constraint1 = 0;
529 GRBLinExpr constraint2 = 0;
531 for (
typename vector<graph_vertex_type>::const_iterator it2 = oddCycle.begin(), it2e = oddCycle.end(); it2 != it2e; ++it2)
533 graph_vertex_type u = *it2;
534 constraint1 += vColorBits[u<<1];
535 constraint2 += vColorBits[(u<<1)+1];
538 sprintf(buf,
"ODD%lu_%u", v, m_constrs_num++);
539 optModel.addConstr(constraint1 >= 1, buf);
540 sprintf(buf,
"ODD%lu_%u", v, m_constrs_num++);
541 optModel.addConstr(constraint1 <= cycleLength-1, buf);
542 sprintf(buf,
"ODD%lu_%u", v, m_constrs_num++);
543 optModel.addConstr(constraint2 >= 1, buf);
544 sprintf(buf,
"ODD%lu_%u", v, m_constrs_num++);
545 optModel.addConstr(constraint2 <= cycleLength-1, buf);
551 template <
typename GraphType>
557 #ifdef DEBUG_LPCOLORING
558 sprintf(buf,
"%u.lp", m_lp_iters);
562 int32_t optStatus = optModel.get(GRB_IntAttr_Status);
563 if (optStatus == GRB_INFEASIBLE)
566 sprintf(buf,
"%u.lp", m_lp_iters);
569 optModel.computeIIS();
570 sprintf(buf,
"%u.ilp", m_lp_iters);
573 assert_msg(optStatus != GRB_INFEASIBLE,
"model is infeasible");
578 template <
typename GraphType>
581 #ifdef DEBUG_LPCOLORING
584 vector<GRBVar> vColorBits;
585 vector<GRBVar> vEdgeBits;
591 env.set(GRB_IntParam_OutputFlag, 0);
594 env.set(GRB_IntParam_Threads, this->m_threads);
596 env.set(GRB_IntParam_Method, -1);
598 GRBModel optModel = GRBModel(env);
601 set_optimize_model(vColorBits, vEdgeBits, obj, optModel);
602 #ifndef DEBUG_NOANCHOR
603 set_anchor(vColorBits);
606 solve_model(optModel);
608 #ifdef DEBUG_LPCOLORING
609 printf(
"\nLP %u solution: ", m_lp_iters);
610 print_solution(vColorBits);
615 non_integer_num(vColorBits, vEdgeBits, curInfo);
616 #ifdef DEBUG_NONINTEGERS
626 adjust_variable_pair_in_objective(vColorBits, obj);
631 optModel.setObjective(obj);
635 add_odd_cycle_constraints(vColorBits, optModel);
637 solve_model(optModel);
639 #ifdef DEBUG_LPCOLORING
640 printf(
"LP %u solution: ", m_lp_iters);
641 print_solution(vColorBits);
645 non_integer_num(vColorBits, vEdgeBits, curInfo);
646 #ifdef DEBUG_NONINTEGERS
655 #ifdef DEBUG_NONINTEGERS
658 vLPNumIter.push_back(m_lp_iters);
664 apply_solution(vColorBits);
668 #ifdef DEBUG_LPCOLORING
671 return this->calc_cost(this->m_vColor);
674 template <
typename GraphType>
677 for (uint32_t i = 0, ie = this->m_vColor.size(); i != ie; ++i)
679 GRBVar
const& var1 = vColorBits[i<<1];
680 GRBVar
const& var2 = vColorBits[(i<<1)+1];
681 int8_t b1 = round(var1.get(GRB_DoubleAttr_X));
682 int8_t b2 = round(var2.get(GRB_DoubleAttr_X));
684 this->m_vColor[i] =
std::min((b1<<1)+b2, (int8_t)this->color_num()-1);
690 template <
typename GraphType>
695 non_integer_num(vColorBits, vEdgeBits, curInfo);
699 bool modifyFlag =
false;
700 for (uint32_t i = 0, ie = vColorBits.size(); i < ie; i += 2)
702 GRBVar
const& var1 = vColorBits[i];
703 GRBVar
const& var2 = vColorBits[i+1];
704 double value1 = var1.get(GRB_DoubleAttr_X);
705 double value2 = var2.get(GRB_DoubleAttr_X);
707 if (!(value1 == 0.5 && value2 == 0.5))
710 GRBColumn column[2] = {
711 optModel.getCol(var1),
712 optModel.getCol(var2)
719 bool mValidColorBits[2][2] = {{
true,
true}, {
true,
true}};
720 if (this->color_num() == base_type::THREE)
721 mValidColorBits[1][1] =
false;
723 uint32_t invalidCount = 0;
724 bool failFlag =
false;
726 for (uint32_t j = 0; j != 2 && !failFlag; ++j)
727 for (uint32_t k = 0, ke = column[j].size(); k != ke; ++k)
729 GRBConstr constr = column[j].getConstr(k);
733 char sense = constr.get(GRB_CharAttr_Sense);
734 curConstrInfo[0].
set(optModel.getCoeff(constr, var1), sense);
735 curConstrInfo[1].
set(optModel.getCoeff(constr, var2), sense);
738 if (!curConstrInfo[0].same_direction(prevConstrInfo[0]) || !curConstrInfo[1].
same_direction(prevConstrInfo[1]))
745 for (int32_t b1 = 0; b1 != 2; ++b1)
746 for (int32_t b2 = 0; b2 != 2; ++b2)
748 if (mValidColorBits[b1][b2])
750 double delta = curConstrInfo[0].
coeff*(b1-value1) + curConstrInfo[1].coeff*(b2-value2);
751 if ((sense ==
'>' && delta < 0.0) || (sense ==
'<' && delta > 0.0))
753 mValidColorBits[b1][b2] =
false;
760 if (invalidCount == 4)
766 prevConstrInfo[0] = curConstrInfo[0];
767 prevConstrInfo[1] = curConstrInfo[1];
773 for (int32_t b1 = 0; b1 != 2; ++b1)
774 for (int32_t b2 = 0; b2 != 2; ++b2)
775 if (mValidColorBits[b1][b2])
777 vColorBits[i].
set(GRB_DoubleAttr_UB, b1);
778 vColorBits[i].set(GRB_DoubleAttr_LB, b1);
779 vColorBits[i+1].set(GRB_DoubleAttr_UB, b2);
780 vColorBits[i+1].set(GRB_DoubleAttr_LB, b2);
787 if (!modifyFlag)
break;
789 solve_model(optModel);
792 non_integer_num(vColorBits, vEdgeBits, curInfo);
797 template <
typename GraphType>
801 if (!this->has_precolored())
805 edge_iterator_type ei, eie;
809 for (boost::tie(ei, eie) = boost::edges(this->m_graph); ei != eie; ++ei)
810 count += refine_color(*ei);
818 template <
typename GraphType>
821 graph_vertex_type v[2] = {
822 boost::source(e, this->m_graph),
823 boost::target(e, this->m_graph)
830 int8_t bestColor[2] = {0, 0};
833 for (c[0] = 0; c[0] != this->color_num(); c[0] += 1)
834 for (c[1] = 0; c[1] != this->color_num(); c[1] += 1)
836 edge_weight_type curCost = 0;
837 typename boost::graph_traits<graph_type>::adjacency_iterator ui, uie;
838 for (int32_t i = 0; i != 2; ++i)
842 graph_vertex_type cv = v[i], ov = v[!i];
843 for (boost::tie(ui, uie) = boost::adjacent_vertices(cv, this->m_graph); ui != uie; ++ui)
845 graph_vertex_type u = *ui;
846 if (u != ov && c[i] == this->m_vColor[u])
848 std::pair<graph_edge_type, bool> eU2cv = boost::edge(cv, u, this->m_graph);
849 #ifdef DEBUG_LPCOLORING
850 assert(eU2cv.second);
853 curCost +=
std::max((edge_weight_type)1, boost::get(boost::edge_weight, this->m_graph, eU2cv.first));
858 curCost +=
std::max((edge_weight_type)1, boost::get(boost::edge_weight, this->m_graph, e));
859 if (curCost < bestCost)
866 bool retFlag =
false;
867 if (this->m_vColor[v[0]] != bestColor[0] || this->m_vColor[v[1]] != bestColor[1])
869 this->m_vColor[v[0]] = bestColor[0];
870 this->m_vColor[v[1]] = bestColor[1];
880 template <
typename GraphType>
886 #ifdef DEBUG_LPCOLORING
892 template <
typename GraphType>
897 for (vector<GRBVar>::const_iterator it = vVariables.begin(), ite = vVariables.end(); it != ite; ++it)
899 double value = it->get(GRB_DoubleAttr_X);
900 if (value != 0.0 && value != 1.0)
910 template <
typename GraphType>
913 uint32_t precoloredNum = 0;
914 for (
typename vector<graph_vertex_type>::const_iterator vi = vVertex.begin(), vie = vVertex.end(); vi != vie; ++vi)
915 if (this->m_vColor[*vi] >= 0)
917 return precoloredNum;
920 template <
typename GraphType>
923 const char* prefix =
"";
924 for (uint32_t i = 0, ie = vColorBits.size(); i < ie; i += 2)
926 printf(
"%sv%u=(%g,%g)", prefix, i>>1, vColorBits[i].
get(GRB_DoubleAttr_X), vColorBits[i+1].
get(GRB_DoubleAttr_X));
hasher class for graph_edge_type
void non_integer_num(vector< GRBVar > const &vColorBits, vector< GRBVar > const &vEdgeBits, NonIntegerInfo &info) const
void print_solution(vector< GRBVar > const &vColorBits) const
boost::dynamic_bitset m_vVertexHandledByOddCycle
members
void set(double c, char s)
void set_anchor(vector< GRBVar > &vColorBits) const
void write_graph(std::ofstream &out, GraphType const &g, VertexLabelType const &vl, EdgeLabelType const &el)
write graph to graphviz format and convert to pdf. Although Boost.Graph has write_graphviz component...
bool refine_color(graph_edge_type const &e)
bool is_integer(string const &s)
check whether string represents an integer
uint32_t vertex_half_integer_num
number of vertices with half-integer solutions
uint32_t m_constrs_num
record number of constraints
ConstrVariableInfo()
constructor
void add_odd_cycle_constraints(vector< GRBVar > const &vColorBits, GRBModel &optModel)
odd cycle trick from Prof. Baldick
bool is_integer(double value) const
uint32_t edge_non_integer_num
number of edges with non-integer solutions
void initialize()
create the NoStitchGraph (this->m_graph) from the (m_conflict_graph)
NonIntegerInfo()
constructor
uint32_t m_lp_iters
record lp iterations
void adjust_variable_pair_in_objective(vector< GRBVar > const &vColorBits, GRBLinExpr &obj) const
tune objective for each color bit pair of vertex
void set_optimize_model(vector< GRBVar > &vColorBits, vector< GRBVar > &vEdgeBits, GRBLinExpr &obj, GRBModel &optModel)
std::iterator_traits< Iterator >::value_type max(Iterator first, Iterator last)
get max of an array
uint32_t post_refinement()
greedy final coloring refinement
void adjust_conflict_edge_vertices_in_objective(vector< GRBVar > const &vColorBits, GRBLinExpr &obj) const
tune objective for each color bit pair along conflict edges
base class for all graph coloring algorithms
uint32_t lp_iters() const
records the information of non-integer values
LPColoring(graph_type const &g)
constructor
void apply_solution(vector< GRBVar > const &vColorBits)
void rounding_with_binding_analysis(GRBModel &optModel, vector< GRBVar > &vColorBits, vector< GRBVar > &vEdgeBits)
Check string is integer, floating point, number... Convert string to upper/lower cases.
uint32_t vertex_non_integer_num
number of vertices with non-integer solutions
void get_odd_cycles(graph_vertex_type const &v, vector< vector< graph_vertex_type > > &vOddCyle)
std::iterator_traits< Iterator >::value_type min(Iterator first, Iterator last)
get min of an array
bool same_direction(ConstrVariableInfo const &rhs) const
graph_vertex_type max_degree_vertex() const
void solve_model(GRBModel &optModel)
solve model
information for a variable of a constraint
#define assert_msg(condition, message)
assertion with message
uint32_t check_precolored_num(vector< LPColoring< GraphType >::graph_vertex_type > const &vVertex) const
uint32_t edge_half_integer_num
number of edges with half-integer solutions