## Plan to implement crossover operator GPX (graph coloring)

In the last few posts, I described my attempts at using Tabu Search only to find a 49-coloring for the standard DIMACS challenge graph: DSJC500.5 or DSJC500.5.col with 500 vertices and a density of 0.5 . For a graph with V vertices, the density is (d^bar)/(V-1) , where d^bar is the average degree of a vertex in the graph. For a complete graph on V vertices, the density is 1.000 . The graph DSJC500.5.col has 500 vertices and 62624 edges. That means the average degree d^bar is 250.496 and that the density is d^bar over 499, or 0.50199599.

My methods to vary Tabu search were:

(a) Start with brand new random approximate solutions, unrelated to previous approximate solutions, as evaluated by the number of edges in which the two vertices are colored the same (also called a “conflict”).

(b) To let the number of iterations in Tabu search (before re-setting the approximate solution) grow: first to 20, then 30, 40, 50, 90 million then finally 1000 billion iterations.

For example, my latest run searching for a 49-coloring (without success):

First select a random 49 coloring of the 500 vertices…

select = 16 num_try = 1

9 conflicts and CN = 20 at 178 iterations

8 conflicts and CN = 17 at 186 iterations

7 conflicts and CN = 16 at 74873 iterations

6 conflicts and CN = 14 at 118666 iterations

5 conflicts and CN = 12 at 585106 iterations

4 conflicts and CN = 10 at 771607 iterations

3 conflicts and CN = 8 at 2800233 iterations

2 conflicts and CN = 6 at 141303044 iterations

The approximate solution was at “2 edges wrong” after about 141 million moves.

I’ve learnt that 2 conflicts isn’t necessarily closer in time to a solution ( = zero conflicts) than “1 conflict” is. Counting the least number of conflicts reached is a “false evaluator” of average time before a real solution (conflict-free) is found.

Published research on heuristics and algorithms for graph coloring includes genetic algorithm methods. In the hybrid method of Galinier and Hao (1999) a tabu search method is combined with a crossover operator GPX (graph partition crossover), the last of which “genetically combines” two aproximate solutions to produce a child approximate solution. The tabu search method carefully chooses one vertex with a conflict, and changes its color to a new one. It does this over and over. This can mean temporarily moving to a “higher point”: an approximate solution with more conflicts. (compare with method of steepest descent in numerical analysis: actually, “Gradient descent” at Wikipedia: https://en.wikipedia.org/wiki/Gradient_descent ).

In these discrete optimization problems, the “landscape” is very irregular and bumpy. Tabu search chooses an optimal one-color switch, but declares reversing the switch (once it had been done to vertex number 56 from blue to red, to be TABOO (forbidden) in the reverse direction, ? to blue on vertex 56, for the next T moves. T can be a constant: 7, 10, 32. It can depend on the number of conflicts. It can depend in part on random numbers. The first recommendation was T = 10. This avoids looping around in a local minimum, or at least that is the aim. Otherwise, as I have seen, one can enter a periodic trajectory (or orbit) that keeps the number of edges with conflicts always at 1 or more.

My work in mathematics and scientific computation has never brought me to implement genetic algorithms. For me, genetic algorithms always had an aura of voodoo or black magic…

Despite knowing that they are effective in graph coloring, what is perhaps most perplexing is: why do they work? Why does adding a crossover operation yield better algorithms?

Nevertheless, I see now that genetic algorithms have a pedigree of a few decades roughly:

https://en.wikipedia.org/wiki/Genetic_algorithm

(recommended for non-initiates, to get a general feel).

The particular method I’m thinking of is HEAD or the simpler HEAD’ of Moalic and Gondran, which incorporates the GPX crossover operator of Galinier and Hao, and works on populations always of size 2, whereas Galinian and Hao recommended populations of about 10 “fittest” solutions (approximate solutions in graph coloring).

Moalic and Gondran:

“Variations on Memetic Algorithms for Graph Coloring Problems” , at arXiv:

https://arxiv.org/abs/1401.2184

====

My interest in algorithms and heuristics for graph coloring stems from the Hadwiger-Nelson problem:

https://en.wikipedia.org/wiki/Hadwiger%E2%80%93Nelson_problem

If one wants to color the Euclidean plane in such a way that points a unit distance apart will always receive distinct colors, then it is known that a 7-coloring known for 50-odd years if not more will do the trick. It is also known that no 3-coloring of plane will do, due to a simple but clever graph now known as the Moser spindle:

https://en.wikipedia.org/wiki/Moser_spindle .

The minimum number of colors to “color the plane” could therefore be 4, 5, 6, or 7.

Earlier in 2017, I did computer experiments to try to find a representation of a given graph in the Euclidean plane where all edges have unit length. The method involves imagining all edges to have springs covering them, which will push apart two points that are closer to each other than the unit distance, and pull together two points linked by an edge when the two points are farther apart than the unit distance.

I described it in the sci.math newsgroup , from January 28 at the Math Forum:

http://mathforum.org/kb/message.jspa?messageID=10074716

The C program attempted, using the adjacency matrix of the 11-vertex Golomb graph, to find a Euclidean plane unit distance representation for this graph, drawn here:

http://mathworld.wolfram.com/GolombGraph.html

Continuing with this general method for 10, 11, 12, 13, 14, and 15 vertices, in time I found a possible/probable 15-vertex, 35-edge unit distance graph . I can only claim “possible/probable”, because the computations, coordinates and distances are all done in finite precision floating point arithmetic. The graph is drawn here, using PostScript :

All graphs found had chromatic number 4, or less.

If one doubles to 30 vertices, the unit distance graph “solver” could still be usable and useful, i.e. it doesn’t appear to take exponential time to get a very good approximation to a unit distance graph representation, although some initial point configurations will converge to a UDG, whereas other initial point configurations will converge to bogus solutions where vertices are at the same point, or in other cases, the distances do not all approach 1 along the edges. But coloring a graph with 30 vertices and 70 edges with the least number of colors would take an immennse amount of time to decide on 2, 3, 4 or however many minimum colors using the brute force method from the Jan 28 code and later versions up to the study of the 15-edge graphs in March 2017 alluded to above.

This motivated me to take up the study of graph coloring algorithms and heuristics.

For the record, I’ll post the source code for the program that allows 10^12 iterations in tabu search for a 49-coloring of the DIMACS challenge DSJC500.5.col graph.

Note: the number of conflicts is down to 1, but I’ve seen “1 conflict” approximate solutions before.

First select a random 49 coloring of the 500 vertices…

select = 16 num_try = 1

9 conflicts and CN = 20 at 178 iterations

8 conflicts and CN = 17 at 186 iterations

7 conflicts and CN = 16 at 74873 iterations

6 conflicts and CN = 14 at 118666 iterations

5 conflicts and CN = 12 at 585106 iterations

4 conflicts and CN = 10 at 771607 iterations

3 conflicts and CN = 8 at 2800233 iterations

2 conflicts and CN = 6 at 141303044 iterations

1 conflicts and CN = 4 at 170724477 iterations

( it’s using 49 colors to try to color the graph DSJC500.5.col ).

Filename = newtabuzzt99a.c

#include <stdio.h>

#include <stdlib.h>

#include <math.h>

#define MAX_VERTEX 1000

#define MAX_COLORS 500

int adj_mat[MAX_VERTEX][MAX_VERTEX];

int Gamma[MAX_VERTEX][MAX_COLORS];

int Fast_Gamma[MAX_VERTEX][MAX_COLORS];

int tabu_list[MAX_VERTEX][MAX_COLORS];

int matches_best_delta[MAX_VERTEX][MAX_COLORS];

static unsigned long long Q[2097152], carry=0;

unsigned long long B64MWC(void)

{ unsigned long long t,x; static int jm=2097151;

jm=(jm+1)&2097151;

x=Q[jm]; t=(x<<28)+carry;

carry=(x>>36)-(t<x);

return (Q[jm]=t-x);

}

#define CNG ( cng=6906969069ULL*cng+13579 )

#define XS ( xs^=(xs<<13), xs^=(xs>>17), xs^=(xs<<43) )

#define KISS ( B64MWC()+CNG+XS )

int main(int argc, char *argv[] )

{

FILE *in;

unsigned long long im,cng=123456789987654321ULL, xs=362436069362436069ULL;

unsigned long long mask = 2147483647ULL;

unsigned long long kiss1, kiss2, seed;

int num_colors;

int num_vertex;

int num_edges;

int num_conflicts;

int k;

int i;

int j;

int sum_adj_mat;

int ru;

int b;

int color;

int init_coloring[MAX_VERTEX];

int old_init_coloring[MAX_VERTEX];

int conflict_count;

int sum_Gamma;

int best_delta;

int delta;

int best_vertex;

int best_color;

long iter_counter;

int A;

float alpha;

int tl;

int b_tabu_rnd;

int tabu_rand;

int color_tabu;

int min_conflicts;

int job_done;

int differences;

FILE *out;

int is_round_one;

int problem_solved;

int cnodes;

int select;

int temp;

long max_iter;

int candidate_count;

int b_candidate;

int can_number;

int num_try;

int ijk;

select = 2;

is_round_one = 1;

out = fopen(“solution407”, “w”);

alpha = 0.60;

A = 10;

job_done = 0;

seed = 10000845109132110031ULL;

seed = 10686169104128798711ULL;

seed = 16082481622012051572ULL;

seed = 18374146278035435545ULL;

seed = 16152130454896673086ULL;

seed = 5797258859902460576ULL;

seed = 9721747250322017459ULL;

num_colors = atoi(argv[2]);

printf(“num_colors = %d\n” , num_colors);

in = fopen(argv[1], “r”);

fscanf(in, “%d %d”, &num_vertex, &num_edges);

printf(“num_vertex = %d\n”, num_vertex);

printf(“num_edges = %d\n”, num_edges);

for(i = 0 ; i < MAX_VERTEX; i++)

{

for(j = 0 ; j < MAX_VERTEX; j++)

{

adj_mat[i][j] = 0;

}

}

for(k = 0; k < num_edges; k++)

{

fscanf(in, “%d %d”, &j, &i);

if( i != j )

{

adj_mat[i-1][j-1] = 1;

adj_mat[j-1][i-1] = 1;

}

if(i == j)

{

printf(“warning: edge %d is = %d %d\n”, k+1, i, j);

}

}

fclose(in);

sum_adj_mat = 0;

for(i = 0 ; i < num_vertex; i++)

{

for(j = 0 ; j < num_vertex; j++)

{

if( (adj_mat[i][j] == 1) && (i!=j) )

{

sum_adj_mat++;

}

}

}

printf(“sum_adj_mat = %d\n”, sum_adj_mat);

printf(“\n”);

for(ijk = 0; ijk < 50; ijk++)

{

job_done = 0;

printf(“ijk = %d\n”, ijk);

num_colors = 500;

num_try = 0;

while( job_done == 0)

{

min_conflicts = 1000000;

b_tabu_rnd = ((int) mask)/A;

b_tabu_rnd = b_tabu_rnd + 1;

printf(“\n”);

printf(“First select a random %d coloring of the %d vertices…\n”, num_colors, num_vertex);

b = ((int) mask)/num_colors;

b = b+1;

if(is_round_one == 1)

{

for(i = 0; i< num_vertex; i++)

{

kiss1 = KISS;

kiss1 = KISS;

kiss1 = KISS;

kiss1 = KISS;

kiss1 = KISS;

kiss1 = KISS;

kiss1 = KISS;

kiss1 = KISS;

kiss1 = KISS;

kiss1 = KISS;

kiss2 = kiss1^seed;

ru = (int) (kiss2&mask);

color = ru/b;

init_coloring[i] = color ;

}

}

if(is_round_one == 0)

{

for(i = 0; i< num_vertex; i++)

{

init_coloring[i] = old_init_coloring[i];

}

select = select + 1;

select = select%num_colors;

for(i = 0; i< num_vertex; i++)

{

if( (init_coloring[i] == select) || (init_coloring[i] == num_colors) )

{

temp = init_coloring[i];

if(temp == select)

{

init_coloring[i]=num_colors;

}

if(temp == num_colors)

{

init_coloring[i]=select ;

}

}

}

for(i = 0; i< num_vertex; i++)

{

if( init_coloring[i] == num_colors )

{

kiss1 = KISS;

kiss2 = kiss1^seed;

ru = (int) (kiss2&mask);

color = ru/b;

init_coloring[i] = color ;

}

}

}

/****************************************************

if(is_round_one == 0)

{

end block

***********************************************/

if(num_colors > 51)

{

max_iter = 100000;

}

if(num_colors == 51)

{

max_iter = 800000;

}

if(num_colors == 50)

{

max_iter = 1000000;

}

if(num_colors <= 49)

{

max_iter = (long) 1000000000000 ;

}

num_try = num_try +1;

for(i = 0 ; i < num_vertex; i++)

{

for(j = 0 ; j < num_colors; j++)

{

tabu_list[i][j] = 0;

}

}

for(iter_counter = 1; iter_counter < max_iter ; iter_counter++)

{

problem_solved = 0;

if( iter_counter == 1)

{

num_conflicts = 0;

printf(“select = %d “, select);

printf(“num_try = %d\n”, num_try);

for(i = 0; i < num_vertex; i++)

{

for(j = 0; j < num_colors; j++)

{

conflict_count = 0;

for(k = 0; k< num_vertex; k++)

{

if( (i != k) && (adj_mat[i][k] == 1) && ( j == init_coloring[k]) )

{

conflict_count++;

}

}

Gamma[i][j] = conflict_count;

Fast_Gamma[i][j] = Gamma[i][j];

}

num_conflicts = num_conflicts + Gamma[i][init_coloring[i]];

}

num_conflicts = num_conflicts/2;

}

if( num_conflicts < min_conflicts )

{

min_conflicts = num_conflicts;

if(num_conflicts < 10)

{

printf(“%d conflicts and CN = %d at %d iterations\n”, num_conflicts, cnodes, iter_counter);

}

if(num_conflicts == 0)

{

problem_solved = 1;

break;

}

}

best_delta = 1000000;

for(i = 0; i < num_vertex; i++)

{

for(j = 0; j < num_colors; j++)

{

matches_best_delta[i][j] = 0;

if(j != init_coloring[i])

{

delta = Fast_Gamma[i][j] – Fast_Gamma[i][init_coloring[i]];

if( delta < best_delta )

{

if( (tabu_list[i][j] < iter_counter) || ((tabu_list[i][j] >= iter_counter)&&((delta + num_conflicts+1) <= min_conflicts )))

{

best_delta = delta;

}

}

}

}

}

candidate_count = 0;

for(i = 0; i < num_vertex; i++)

{

for(j = 0; j < num_colors; j++)

{

if(j != init_coloring[i])

{

delta = Fast_Gamma[i][j] – Fast_Gamma[i][init_coloring[i]];

if( delta == best_delta )

{

if( (tabu_list[i][j] < iter_counter) || ((tabu_list[i][j] >= iter_counter)&&((delta + num_conflicts+1) <= min_conflicts )))

{

matches_best_delta[i][j] = 1;

candidate_count++;

}

}

}

}

}

b_candidate = ((int) mask)/candidate_count;

b_candidate = b_candidate + 1;

ru = (int) (KISS&mask);

can_number = ru/b_candidate;

candidate_count = 0;

for(i = 0; i < num_vertex; i++)

{

for(j = 0; j < num_colors; j++)

{

if(matches_best_delta[i][j] == 1)

{

if(candidate_count == can_number)

{

best_vertex = i;

best_color = j;

}

candidate_count++;

}

}

}

num_conflicts = num_conflicts + best_delta;

cnodes = 0;

for(i = 0; i < num_vertex; i++)

{

if( Fast_Gamma[i][init_coloring[i]] > 0)

{

cnodes++;

}

}

tl = (int) (floor( alpha*((float) cnodes ) + 0.5));

kiss1 = KISS;

kiss2 = kiss1^seed;

ru = (int) (kiss2&mask);

tabu_rand = ru/b_tabu_rnd;

tl = tl + tabu_rand;

tl = 30;

color_tabu = init_coloring[best_vertex];

tabu_list[best_vertex][color_tabu] = tl + iter_counter;

init_coloring[best_vertex] = best_color;

for(i = 0; i < num_vertex; i++)

{

if( adj_mat[i][best_vertex] == 1)

{

Fast_Gamma[i][color_tabu] = Fast_Gamma[i][color_tabu]-1;

Fast_Gamma[i][best_color] = Fast_Gamma[i][best_color]+1;

}

}

}

/*******************************************

for(iter_counter = 1; iter_counter < max_iter ; iter_counter++)

{

end loop …

**************************************/

if( (problem_solved == 1) )

{

if(num_colors <= 50)

{

fprintf(out, “solution for %d colors:\n”, num_colors);

for(i = 0; i < num_vertex; i++ )

{

fprintf(out, “%d %d\n”, i+1, init_coloring[i]+1);

}

fprintf(out, “\n”);

fflush(out);

}

num_colors = num_colors – 1;

// printf(“num_colors = %d\n”, num_colors);

num_try = 0;

is_round_one = 0;

for(i = 0; i < num_vertex; i++ )

{

old_init_coloring[i] = init_coloring[i];

}

}

if( (problem_solved == 0)&& (num_colors == 50)&&(num_try >=1) )

{

job_done = 1;

}

if( (problem_solved == 0)&& (num_colors == 49)&&(num_try >=1) )

{

job_done = 1;

}

}

/*******************************************************

while( job_done == 0)

{

end loop …

********************************************************/

}

/********************************************

for(ijk = 0; ijk < 50; ijk++)

{

end loop …

*****************************************/

fclose(out);

return 0;

}