DSDP
color.c
Go to the documentation of this file.
1 #include "dsdp5.h"
2 #include <string.h>
3 #include <math.h>
4 #include <stdio.h>
5 #include <stdlib.h>
6 
12 char help2[]="\nA positive semidefinite relaxation of the\n\
13 graph k coloring problem can be rewritten as\n\n\
14  Find X>=0 \n\
15  such that X_ij <= 1 - 1/(k-1) for all edges (i,j).\n\
16 ";
17 
18 char help[]="DSDP Usage: color <graph filename> ";
19 
20 static int ReadGraph(char*,int *, int *,int**, int **, double **);
21 static int ParseGraphline(char*,int*,int*,double*,int*);
22 static int RandomizedColor(DSDP, SDPCone, int, int[], int[], int);
23 int MinColoring(int argc,char *argv[]);
24 
25 
26 int main(int argc,char *argv[]){
27  int info;
28  info=MinColoring(argc,argv);
29  return 0;
30 }
31 
39 int MinColoring(int argc,char *argv[]){
40 
41  int i,kk,vari,info;
42  int *node1,*node2,nedges,nnodes;
43  int *iptr1,*iptr2;
44  double *weight,*yy1,*yy2,bb;
45  DSDPTerminationReason reason;
46  SDPCone sdpcone;
47  BCone bcone;
48  DSDP dsdp;
49 
50  if (argc<2){ printf("%s\n%s",help2,help); return(1); }
51 
52  info = ReadGraph(argv[1],&nnodes,&nedges,&node1,&node2,&weight);
53  if (info){ printf("Problem reading file\n"); return 1; }
54 
55  info = DSDPCreate(nnodes+nedges,&dsdp);
56 
57  info = DSDPCreateSDPCone(dsdp,1,&sdpcone);
58  info = SDPConeSetBlockSize(sdpcone,0,nnodes);
59  info = SDPConeSetSparsity(sdpcone,0,nnodes+nedges+1);
60 
61  info = DSDPCreateBCone(dsdp,&bcone);
62 
63  if (info){ printf("Out of memory\n"); return 1; }
64 
65 
66  /* Formulate the problem from the data */
67  /* Create data matrices */
68 
69  /* Diagonal elements of X(i,i) must equal 1.0 */
70  iptr1=(int*)malloc(nnodes*sizeof(int));
71  yy1=(double*)malloc(nnodes*sizeof(double));
72  for (i=0;i<nnodes;i++){
73  iptr1[i]=(i+2)*(i+1)/2-1;
74  yy1[i]=1.0;
75  }
76  for (i=0;i<nnodes;i++){
77  info=SDPConeSetASparseVecMat(sdpcone,0,i+1,nnodes,1.0,0,iptr1+i,yy1+i,1);
78  if (info) printf("ERROR 1: %d \n",i);
79  info=DSDPSetDualObjective(dsdp,i+1,1.0);
80  }
81 
82  /* For each nonzero element (i,j) of the matrix, X(i,j) must be less than 1 - 1/nnodes */
83  bb=2-2.0/nnodes;
84  iptr2=(int*)malloc(nedges*sizeof(int));
85  yy2=(double*)malloc(nedges*sizeof(double));
86  for (i=0;i<nedges;i++){
87  iptr2[i]=(node1[i])*(node1[i]+1)/2+node2[i];
88  yy2[i]=1.0;
89  }
90  info = BConeAllocateBounds(bcone,nedges);
91  for (i=0; i<nedges; i++){
92  vari=nnodes+i+1;
93  info = SDPConeSetSparseVecMat(sdpcone,0,vari,nnodes,0,iptr2+i,yy2+i,1);
94  if (info) printf("ERROR 2: %d %d \n",i,vari);
95  info = BConeSetPSlackVariable(bcone,vari);
96  if (info) printf("ERROR 3: %d %d \n",i,vari);
97  info = DSDPSetDualObjective(dsdp,vari,bb);
98  }
99 
100 
101  /* Get read to go */
102  info=DSDPSetPotentialParameter(dsdp,5);
103 
104  for (kk=1; kk<argc-1; kk++){
105  if (strncmp(argv[kk],"-dloginfo",8)==0){
106  info=DSDPLogInfoAllow(DSDP_TRUE,0);
107  } else if (strncmp(argv[kk],"-params",7)==0){
108  info=DSDPReadOptions(dsdp,argv[kk+1]);
109  } else if (strncmp(argv[kk],"-help",7)==0){
110  printf("%s\n",help);
111  }
112  }
113  info=DSDPSetOptions(dsdp,argv,argc);
114 
115  if (info){ printf("Out of memory\n"); return 1; }
116  info = DSDPSetStandardMonitor(dsdp,1);
117 
118  info = DSDPSetup(dsdp);
119  if (info){ printf("Out of memory\n"); return 1; }
120 
121  info = DSDPSolve(dsdp);
122  if (info){ printf("Numerical error\n"); return 1; }
123  info = DSDPStopReason(dsdp,&reason);
124 
125  if (reason!=DSDP_INFEASIBLE_START){ /* Randomized solution strategy */
126  info=RandomizedColor(dsdp, sdpcone, nnodes, node1, node2, nedges);
127  }
128 
129  info = DSDPDestroy(dsdp);
130 
131  free(node1);free(node2);free(weight);
132  free(iptr1);
133  free(yy1);
134  free(iptr2);
135  free(yy2);
136 
137  return 0;
138 }
139 
140 static int GetXRow(double xmat[],double xrow[],int row,int n){
141  int i,i1=row*(row+1)/2;
142  for (i=0;i<row;i++){xrow[i]=xmat[i1+i];}
143  for (i=row;i<n;i++){xrow[i]=xmat[i1+row];i1+=i+1;}
144  return 0;
145 }
146 
147 typedef struct {
148  int index;double val;
149 } orderVec;
150 
151 static int cut_comp(const void *e1,const void *e2){ /* Used in qsort routine */
152  double d1=((orderVec*)e1)->val, d2=((orderVec*)e2)->val;
153  if (d1<d2) return (1);
154  else if (d1>d2) return (-1);
155  return(0);
156 }
157 
158 static int RemoveNode(int node, int node1[], int node2[], int *nedges){
159  int i,nnedges=*nedges;
160  for (i=0;i<nnedges;i++){
161  if (node1[i]==node || node2[i]==node){
162  node1[i]=node1[nnedges-1];
163  node2[i]=node2[nnedges-1];
164  nnedges--;
165  if (i < nnedges) i--;
166  }
167  }
168  *nedges=nnedges;
169  return 0;
170 }
171 
172 static int Connected(int n1, int n2, int node1[], int node2[], int nedges){
173  int i;
174  if (n1==n2) return 1;
175  for (i=0;i<nedges;i++){
176  if (node1[i]==n1 && node2[i]==n2){ return 1;}
177  if (node1[i]==n2 && node2[i]==n1){ return 1;}
178  }
179  return 0;
180 }
181 
182 static int HighDegree(int node1[], int node2[], int nedges, int iwork[], int nnodes){
183  int i,nmax=0,maxdegree=-1;
184  for (i=0;i<nnodes;i++) iwork[i]=0;
185  for (i=0;i<nedges;i++){
186  iwork[node1[i]]++; iwork[node2[i]]++;
187  }
188  for (i=0;i<nnodes;i++){ if (iwork[i]>maxdegree){nmax=i; maxdegree=iwork[i];} }
189  return nmax;
190 }
191 
192 static int First(int coloring[], int nnodes){
193  int i,nmax=nnodes;
194  for (i=0;i<nnodes;i++){
195  if (coloring[i]==0){
196  nmax=i; return nmax;
197  }
198  }
199  return -1;
200 }
201 
202 static int RandomizedColor(DSDP dsdp, SDPCone sdpcone, int nnodes, int node1[], int node2[], int nedges){
203  int i,j,nodek,nn,info,flag,coloring=0,maxvertex;
204  int *degree,*color,*cgroup,ngsize,uncolored=nnodes;
205  int tnedges=nedges;
206  double *xrow,*xptr;
207  orderVec *vorder;
208 
209  xrow=(double*)malloc(nnodes*sizeof(double));
210  color=(int*)malloc(nnodes*sizeof(int));
211  cgroup=(int*)malloc(nnodes*sizeof(int));
212  degree=(int*)malloc(nnodes*sizeof(int));
213  vorder=(orderVec*)malloc(nnodes*sizeof(orderVec));
214 
215  for (i=0;i<nnodes;i++){ color[i]=0;}
216  info=DSDPComputeX(dsdp);
217  info=SDPConeGetXArray(sdpcone,0,&xptr,&nn);
218 
219  while (uncolored>0){
220 
221  coloring++;
222 
223  maxvertex=First(color,nnodes);
224  maxvertex=HighDegree(node1,node2,tnedges,degree,nnodes);
225 
226  cgroup[0]=maxvertex;ngsize=1;
227 
228  info=GetXRow(xptr,xrow,maxvertex,nnodes);
229 
230  for (i=0;i<nnodes;i++){vorder[i].index=i; vorder[i].val = xrow[i];}
231  qsort( (void*)vorder, nnodes, sizeof(orderVec), cut_comp);
232 
233  for (i=0;i<nnodes;i++){
234  nodek=vorder[i].index;
235  if (color[nodek]==0){
236  for (flag=0,j=0;j<ngsize;j++){
237  if (Connected(nodek,cgroup[j],node1,node2,tnedges) ){flag=1;}
238  }
239  if (flag==0){ cgroup[ngsize]=nodek; ngsize++; }
240  }
241  }
242  for (i=0;i<ngsize;i++){
243  color[cgroup[i]]=coloring; uncolored--;
244  RemoveNode(cgroup[i],node1,node2,&tnedges);
245  }
246  }
247  printf("\nCOLORS: %d\n",coloring);
248  free(xrow);
249  free(color);
250  free(cgroup);
251  free(degree);
252  free(vorder);
253  return(0);
254 }
255 
256 
257 #define BUFFERSIZ 100
258 
259 #undef __FUNCT__
260 #define __FUNCT__ "ParseGraphline"
261 static int ParseGraphline(char * thisline,int *row,int *col,double *value,
262  int *gotem){
263 
264  int temp;
265  int rtmp,coltmp;
266 
267  rtmp=-1, coltmp=-1, *value=0.0;
268  temp=sscanf(thisline,"%d %d %lf",&rtmp,&coltmp,value);
269  if (temp==3 && rtmp>0 && coltmp>0) *gotem=3;
270  else if (temp==2 && rtmp>0 && coltmp>0){ *value = 1.0; *gotem=3;}
271  else *gotem=0;
272  *row=rtmp-1; *col=coltmp-1;
273 
274  return(0);
275 }
276 
277 
278 #undef __FUNCT__
279 #define __FUNCT__ "ReadGraph"
280 int ReadGraph(char* filename,int *nnodes, int *nedges,
281  int**n1, int ** n2, double **wght){
282 
283  FILE*fp;
284  char thisline[BUFFERSIZ]="*";
285  int i,k=0,line=0,nodes,edges,gotone=3;
286  int *node1,*node2;
287  double *weight;
288  int info,row,col;
289  double value;
290 
291  fp=fopen(filename,"r");
292  if (!fp){printf("Cannot open file %s !",filename);return(1);}
293 
294  while(!feof(fp) && (thisline[0] == '*' || thisline[0] == '"') ){
295  fgets(thisline,BUFFERSIZ,fp); line++;
296  }
297 
298  if (sscanf(thisline,"%d %d",&nodes, &edges)!=2){
299  printf("First line must contain the number of nodes and number of edges\n");
300  return 1;
301  }
302 
303  node1=(int*)malloc(edges*sizeof(int));
304  node2=(int*)malloc(edges*sizeof(int));
305  weight=(double*)malloc(edges*sizeof(double));
306 
307  for (i=0; i<edges; i++){ node1[i]=0;node2[i]=0;weight[i]=0.0;}
308 
309  while(!feof(fp) && gotone){
310  thisline[0]='\0';
311  fgets(thisline,BUFFERSIZ,fp); line++;
312  info = ParseGraphline(thisline,&row,&col,&value,&gotone);
313  if (gotone && value!=0.0 && k<edges &&
314  col < nodes && row < nodes && col >= 0 && row >= 0){
315  if (row<col){info=row;row=col;col=info;}
316  if (row == col){}
317  else {
318  node1[k]=row; node2[k]=col;
319  weight[k]=value; k++;
320  }
321  } else if (gotone && k>=edges) {
322  printf("To many edges in file.\nLine %d, %s\n",line,thisline);
323  return 1;
324  } else if (gotone&&(col >= nodes || row >= nodes || col < 0 || row < 0)){
325  printf("Invalid node number.\nLine %d, %s\n",line,thisline);
326  return 1;
327  }
328  }
329  *nnodes=nodes; *nedges=edges;
330  *n1=node1; *n2=node2; *wght=weight;
331  return 0;
332 }
int DSDPCreate(int m, DSDP *dsdpnew)
Create a DSDP solver. FIRST DSDP routine!
Definition: dsdpsetup.c:30
int DSDPDestroy(DSDP dsdp)
Free the internal data structures of the solver and the cones associated with it.
Definition: dsdpsetup.c:496
int DSDPReadOptions(DSDP dsdp, char filename[])
Read DSDP parameters from a file.
int DSDPCreateBCone(DSDP dsdp, BCone *dspcone)
Create a new cone that represents bounds on the y variables.
Definition: dbounds.c:467
int DSDPSetOptions(DSDP dsdp, char *runargs[], int nargs)
Read command line arguments to set options in DSDP.
int DSDPSetDualObjective(DSDP dsdp, int i, double bi)
Set the objective vector b in (D).
Definition: dsdpsetdata.c:25
Internal structures for the DSDP solver.
Definition: dsdp.h:65
DSDPTerminationReason
There are many reasons to terminate the solver.
int SDPConeGetXArray(SDPCone sdpcone, int blockj, double *xx[], int *nn)
After applying the solver, set a pointer to the array in the object with the solution X.
Definition: dsdpadddata.c:328
The API to DSDP for those applications using DSDP as a subroutine library.
int DSDPSetup(DSDP dsdp)
Set up data structures in the solver and the cones associated with it.
Definition: dsdpsetup.c:193
int BConeAllocateBounds(BCone bcone, int nnz)
Set a surplus variable in constraint in (P).
Definition: dbounds.c:645
int BConeSetPSlackVariable(BCone bcone, int vari)
Set a slack variable to a constraint in (P).
Definition: dbounds.c:607
int DSDPSetPotentialParameter(DSDP dsdp, double rho)
Set the potential parameter.
Definition: dsdpsetdata.c:765
int DSDPStopReason(DSDP dsdp, DSDPTerminationReason *reason)
Copy the reason why the solver terminated.
Definition: dsdpsetdata.c:582
int MinColoring(int argc, char *argv[])
SDP relaxation of k-coloring problem.
Definition: color.c:39
int DSDPSolve(DSDP dsdp)
Apply DSDP to the problem.
Definition: dsdpsetup.c:343
int SDPConeSetBlockSize(SDPCone sdpcone, int blockj, int n)
Set the dimension of one block in the semidefinite cone.
Definition: dsdpadddata.c:535
Internal structure for semidefinite cone.
Definition: dsdpsdp.h:80
int SDPConeSetSparsity(SDPCone sdpcone, int blockj, int nnz)
Set the number of nonzero matrices in a block of the semidefinite cone.
Definition: dsdpadddata.c:596
int DSDPComputeX(DSDP dsdp)
Compute the X variables.
Definition: dsdpx.c:55
int SDPConeSetASparseVecMat(SDPCone sdpcone, int blockj, int vari, int n, double alpha, int ishift, const int ind[], const double val[], int nnz)
Set data matrix in a sparse format.
struct BCone_C * BCone
The BCone object points to lower and upper bounds on the variable y in (D).
Definition: dsdp5.h:28
int DSDPSetStandardMonitor(DSDP dsdp, int k)
Print at every kth iteration.
Definition: dsdpprintout.c:153