You can not select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
242 lines
7.3 KiB
242 lines
7.3 KiB
/* wclique.c (maximum weight clique, Ostergard's algorithm) */
|
|
|
|
/***********************************************************************
|
|
* This code is part of GLPK (GNU Linear Programming Kit).
|
|
*
|
|
* Two subroutines sub() and wclique() below are intended to find a
|
|
* maximum weight clique in a given undirected graph. These subroutines
|
|
* are slightly modified version of the program WCLIQUE developed by
|
|
* Patric Ostergard <http://www.tcs.hut.fi/~pat/wclique.html> and based
|
|
* on ideas from the article "P. R. J. Ostergard, A new algorithm for
|
|
* the maximum-weight clique problem, submitted for publication", which
|
|
* in turn is a generalization of the algorithm for unweighted graphs
|
|
* presented in "P. R. J. Ostergard, A fast algorithm for the maximum
|
|
* clique problem, submitted for publication".
|
|
*
|
|
* USED WITH PERMISSION OF THE AUTHOR OF THE ORIGINAL CODE.
|
|
*
|
|
* Changes were made by Andrew Makhorin <mao@gnu.org>.
|
|
*
|
|
* GLPK is free software: you can redistribute it and/or modify it
|
|
* under the terms of the GNU General Public License as published by
|
|
* the Free Software Foundation, either version 3 of the License, or
|
|
* (at your option) any later version.
|
|
*
|
|
* GLPK is distributed in the hope that it will be useful, but WITHOUT
|
|
* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
|
|
* or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
|
|
* License for more details.
|
|
*
|
|
* You should have received a copy of the GNU General Public License
|
|
* along with GLPK. If not, see <http://www.gnu.org/licenses/>.
|
|
***********************************************************************/
|
|
|
|
#include "env.h"
|
|
#include "wclique.h"
|
|
|
|
/***********************************************************************
|
|
* NAME
|
|
*
|
|
* wclique - find maximum weight clique with Ostergard's algorithm
|
|
*
|
|
* SYNOPSIS
|
|
*
|
|
* #include "wclique.h"
|
|
* int wclique(int n, const int w[], const unsigned char a[],
|
|
* int ind[]);
|
|
*
|
|
* DESCRIPTION
|
|
*
|
|
* The routine wclique finds a maximum weight clique in an undirected
|
|
* graph with Ostergard's algorithm.
|
|
*
|
|
* INPUT PARAMETERS
|
|
*
|
|
* n is the number of vertices, n > 0.
|
|
*
|
|
* w[i], i = 1,...,n, is a weight of vertex i.
|
|
*
|
|
* a[*] is the strict (without main diagonal) lower triangle of the
|
|
* graph adjacency matrix in packed format.
|
|
*
|
|
* OUTPUT PARAMETER
|
|
*
|
|
* ind[k], k = 1,...,size, is the number of a vertex included in the
|
|
* clique found, 1 <= ind[k] <= n, where size is the number of vertices
|
|
* in the clique returned on exit.
|
|
*
|
|
* RETURNS
|
|
*
|
|
* The routine returns the clique size, i.e. the number of vertices in
|
|
* the clique. */
|
|
|
|
struct csa
|
|
{ /* common storage area */
|
|
int n;
|
|
/* number of vertices */
|
|
const int *wt; /* int wt[0:n-1]; */
|
|
/* weights */
|
|
const unsigned char *a;
|
|
/* adjacency matrix (packed lower triangle without main diag.) */
|
|
int record;
|
|
/* weight of best clique */
|
|
int rec_level;
|
|
/* number of vertices in best clique */
|
|
int *rec; /* int rec[0:n-1]; */
|
|
/* best clique so far */
|
|
int *clique; /* int clique[0:n-1]; */
|
|
/* table for pruning */
|
|
int *set; /* int set[0:n-1]; */
|
|
/* current clique */
|
|
};
|
|
|
|
#define n (csa->n)
|
|
#define wt (csa->wt)
|
|
#define a (csa->a)
|
|
#define record (csa->record)
|
|
#define rec_level (csa->rec_level)
|
|
#define rec (csa->rec)
|
|
#define clique (csa->clique)
|
|
#define set (csa->set)
|
|
|
|
#if 0
|
|
static int is_edge(struct csa *csa, int i, int j)
|
|
{ /* if there is arc (i,j), the routine returns true; otherwise
|
|
* false; 0 <= i, j < n */
|
|
int k;
|
|
xassert(0 <= i && i < n);
|
|
xassert(0 <= j && j < n);
|
|
if (i == j) return 0;
|
|
if (i < j) k = i, i = j, j = k;
|
|
k = (i * (i - 1)) / 2 + j;
|
|
return a[k / CHAR_BIT] &
|
|
(unsigned char)(1 << ((CHAR_BIT - 1) - k % CHAR_BIT));
|
|
}
|
|
#else
|
|
#define is_edge(csa, i, j) ((i) == (j) ? 0 : \
|
|
(i) > (j) ? is_edge1(i, j) : is_edge1(j, i))
|
|
#define is_edge1(i, j) is_edge2(((i) * ((i) - 1)) / 2 + (j))
|
|
#define is_edge2(k) (a[(k) / CHAR_BIT] & \
|
|
(unsigned char)(1 << ((CHAR_BIT - 1) - (k) % CHAR_BIT)))
|
|
#endif
|
|
|
|
static void sub(struct csa *csa, int ct, int table[], int level,
|
|
int weight, int l_weight)
|
|
{ int i, j, k, curr_weight, left_weight, *p1, *p2, *newtable;
|
|
newtable = xcalloc(n, sizeof(int));
|
|
if (ct <= 0)
|
|
{ /* 0 or 1 elements left; include these */
|
|
if (ct == 0)
|
|
{ set[level++] = table[0];
|
|
weight += l_weight;
|
|
}
|
|
if (weight > record)
|
|
{ record = weight;
|
|
rec_level = level;
|
|
for (i = 0; i < level; i++) rec[i] = set[i];
|
|
}
|
|
goto done;
|
|
}
|
|
for (i = ct; i >= 0; i--)
|
|
{ if ((level == 0) && (i < ct)) goto done;
|
|
k = table[i];
|
|
if ((level > 0) && (clique[k] <= (record - weight)))
|
|
goto done; /* prune */
|
|
set[level] = k;
|
|
curr_weight = weight + wt[k];
|
|
l_weight -= wt[k];
|
|
if (l_weight <= (record - curr_weight))
|
|
goto done; /* prune */
|
|
p1 = newtable;
|
|
p2 = table;
|
|
left_weight = 0;
|
|
while (p2 < table + i)
|
|
{ j = *p2++;
|
|
if (is_edge(csa, j, k))
|
|
{ *p1++ = j;
|
|
left_weight += wt[j];
|
|
}
|
|
}
|
|
if (left_weight <= (record - curr_weight)) continue;
|
|
sub(csa, p1 - newtable - 1, newtable, level + 1, curr_weight,
|
|
left_weight);
|
|
}
|
|
done: xfree(newtable);
|
|
return;
|
|
}
|
|
|
|
int wclique(int n_, const int w[], const unsigned char a_[], int ind[])
|
|
{ struct csa csa_, *csa = &csa_;
|
|
int i, j, p, max_wt, max_nwt, wth, *used, *nwt, *pos;
|
|
double timer;
|
|
n = n_;
|
|
xassert(n > 0);
|
|
wt = &w[1];
|
|
a = a_;
|
|
record = 0;
|
|
rec_level = 0;
|
|
rec = &ind[1];
|
|
clique = xcalloc(n, sizeof(int));
|
|
set = xcalloc(n, sizeof(int));
|
|
used = xcalloc(n, sizeof(int));
|
|
nwt = xcalloc(n, sizeof(int));
|
|
pos = xcalloc(n, sizeof(int));
|
|
/* start timer */
|
|
timer = xtime();
|
|
/* order vertices */
|
|
for (i = 0; i < n; i++)
|
|
{ nwt[i] = 0;
|
|
for (j = 0; j < n; j++)
|
|
if (is_edge(csa, i, j)) nwt[i] += wt[j];
|
|
}
|
|
for (i = 0; i < n; i++)
|
|
used[i] = 0;
|
|
for (i = n-1; i >= 0; i--)
|
|
{ max_wt = -1;
|
|
max_nwt = -1;
|
|
for (j = 0; j < n; j++)
|
|
{ if ((!used[j]) && ((wt[j] > max_wt) || (wt[j] == max_wt
|
|
&& nwt[j] > max_nwt)))
|
|
{ max_wt = wt[j];
|
|
max_nwt = nwt[j];
|
|
p = j;
|
|
}
|
|
}
|
|
pos[i] = p;
|
|
used[p] = 1;
|
|
for (j = 0; j < n; j++)
|
|
if ((!used[j]) && (j != p) && (is_edge(csa, p, j)))
|
|
nwt[j] -= wt[p];
|
|
}
|
|
/* main routine */
|
|
wth = 0;
|
|
for (i = 0; i < n; i++)
|
|
{ wth += wt[pos[i]];
|
|
sub(csa, i, pos, 0, 0, wth);
|
|
clique[pos[i]] = record;
|
|
if (xdifftime(xtime(), timer) >= 5.0 - 0.001)
|
|
{ /* print current record and reset timer */
|
|
xprintf("level = %d (%d); best = %d\n", i+1, n, record);
|
|
timer = xtime();
|
|
}
|
|
}
|
|
xfree(clique);
|
|
xfree(set);
|
|
xfree(used);
|
|
xfree(nwt);
|
|
xfree(pos);
|
|
/* return the solution found */
|
|
for (i = 1; i <= rec_level; i++) ind[i]++;
|
|
return rec_level;
|
|
}
|
|
|
|
#undef n
|
|
#undef wt
|
|
#undef a
|
|
#undef record
|
|
#undef rec_level
|
|
#undef rec
|
|
#undef clique
|
|
#undef set
|
|
|
|
/* eof */
|