00001 #include "dsdp5.h"
00002 #include <string.h>
00003 #include <math.h>
00004 #include <stdio.h>
00005 #include <stdlib.h>
00006
00012 char help2[]="\nA positive semidefinite relaxation of the\n\
00013 graph k coloring problem can be rewritten as\n\n\
00014 Find X>=0 \n\
00015 such that X_ij <= 1 - 1/(k-1) for all edges (i,j).\n\
00016 ";
00017
00018 char help[]="DSDP Usage: color <graph filename> ";
00019
00020 static int ReadGraph(char*,int *, int *,int**, int **, double **);
00021 static int ParseGraphline(char*,int*,int*,double*,int*);
00022 static int RandomizedColor(DSDP, SDPCone, int, int[], int[], int);
00023 int MinColoring(int argc,char *argv[]);
00024
00025
00026 int main(int argc,char *argv[]){
00027 int info;
00028 info=MinColoring(argc,argv);
00029 return 0;
00030 }
00031
00039 int MinColoring(int argc,char *argv[]){
00040
00041 int i,kk,vari,info;
00042 int *node1,*node2,nedges,nnodes;
00043 int *iptr1,*iptr2;
00044 double *weight,*yy1,*yy2,bb;
00045 DSDPTerminationReason reason;
00046 SDPCone sdpcone;
00047 BCone bcone;
00048 DSDP dsdp;
00049
00050 if (argc<2){ printf("%s\n%s",help2,help); return(1); }
00051
00052 info = ReadGraph(argv[1],&nnodes,&nedges,&node1,&node2,&weight);
00053 if (info){ printf("Problem reading file\n"); return 1; }
00054
00055 info = DSDPCreate(nnodes+nedges,&dsdp);
00056
00057 info = DSDPCreateSDPCone(dsdp,1,&sdpcone);
00058 info = SDPConeSetBlockSize(sdpcone,0,nnodes);
00059 info = SDPConeSetSparsity(sdpcone,0,nnodes+nedges+1);
00060
00061 info = DSDPCreateBCone(dsdp,&bcone);
00062
00063 if (info){ printf("Out of memory\n"); return 1; }
00064
00065
00066
00067
00068
00069
00070 iptr1=(int*)malloc(nnodes*sizeof(int));
00071 yy1=(double*)malloc(nnodes*sizeof(double));
00072 for (i=0;i<nnodes;i++){
00073 iptr1[i]=(i+2)*(i+1)/2-1;
00074 yy1[i]=1.0;
00075 }
00076 for (i=0;i<nnodes;i++){
00077 info=SDPConeSetASparseVecMat(sdpcone,0,i+1,nnodes,1.0,0,iptr1+i,yy1+i,1);
00078 if (info) printf("ERROR 1: %d \n",i);
00079 info=DSDPSetDualObjective(dsdp,i+1,1.0);
00080 }
00081
00082
00083 bb=2-2.0/nnodes;
00084 iptr2=(int*)malloc(nedges*sizeof(int));
00085 yy2=(double*)malloc(nedges*sizeof(double));
00086 for (i=0;i<nedges;i++){
00087 iptr2[i]=(node1[i])*(node1[i]+1)/2+node2[i];
00088 yy2[i]=1.0;
00089 }
00090 info = BConeAllocateBounds(bcone,nedges);
00091 for (i=0; i<nedges; i++){
00092 vari=nnodes+i+1;
00093 info = SDPConeSetSparseVecMat(sdpcone,0,vari,nnodes,0,iptr2+i,yy2+i,1);
00094 if (info) printf("ERROR 2: %d %d \n",i,vari);
00095 info = BConeSetPSlackVariable(bcone,vari);
00096 if (info) printf("ERROR 3: %d %d \n",i,vari);
00097 info = DSDPSetDualObjective(dsdp,vari,bb);
00098 }
00099
00100
00101
00102 info=DSDPSetPotentialParameter(dsdp,5);
00103
00104 for (kk=1; kk<argc-1; kk++){
00105 if (strncmp(argv[kk],"-dloginfo",8)==0){
00106 info=DSDPLogInfoAllow(DSDP_TRUE,0);
00107 } else if (strncmp(argv[kk],"-params",7)==0){
00108 info=DSDPReadOptions(dsdp,argv[kk+1]);
00109 } else if (strncmp(argv[kk],"-help",7)==0){
00110 printf("%s\n",help);
00111 }
00112 }
00113 info=DSDPSetOptions(dsdp,argv,argc);
00114
00115 if (info){ printf("Out of memory\n"); return 1; }
00116 info = DSDPSetStandardMonitor(dsdp,1);
00117
00118 info = DSDPSetup(dsdp);
00119 if (info){ printf("Out of memory\n"); return 1; }
00120
00121 info = DSDPSolve(dsdp);
00122 if (info){ printf("Numerical error\n"); return 1; }
00123 info = DSDPStopReason(dsdp,&reason);
00124
00125 if (reason!=DSDP_INFEASIBLE_START){
00126 info=RandomizedColor(dsdp, sdpcone, nnodes, node1, node2, nedges);
00127 }
00128
00129 info = DSDPDestroy(dsdp);
00130
00131 free(node1);free(node2);free(weight);
00132 free(iptr1);
00133 free(yy1);
00134 free(iptr2);
00135 free(yy2);
00136
00137 return 0;
00138 }
00139
00140 static int GetXRow(double xmat[],double xrow[],int row,int n){
00141 int i,i1=row*(row+1)/2;
00142 for (i=0;i<row;i++){xrow[i]=xmat[i1+i];}
00143 for (i=row;i<n;i++){xrow[i]=xmat[i1+row];i1+=i+1;}
00144 return 0;
00145 }
00146
00147 typedef struct {
00148 int index;double val;
00149 } orderVec;
00150
00151 static int cut_comp(const void *e1,const void *e2){
00152 double d1=((orderVec*)e1)->val, d2=((orderVec*)e2)->val;
00153 if (d1<d2) return (1);
00154 else if (d1>d2) return (-1);
00155 return(0);
00156 }
00157
00158 static int RemoveNode(int node, int node1[], int node2[], int *nedges){
00159 int i,nnedges=*nedges;
00160 for (i=0;i<nnedges;i++){
00161 if (node1[i]==node || node2[i]==node){
00162 node1[i]=node1[nnedges-1];
00163 node2[i]=node2[nnedges-1];
00164 nnedges--;
00165 if (i < nnedges) i--;
00166 }
00167 }
00168 *nedges=nnedges;
00169 return 0;
00170 }
00171
00172 static int Connected(int n1, int n2, int node1[], int node2[], int nedges){
00173 int i;
00174 if (n1==n2) return 1;
00175 for (i=0;i<nedges;i++){
00176 if (node1[i]==n1 && node2[i]==n2){ return 1;}
00177 if (node1[i]==n2 && node2[i]==n1){ return 1;}
00178 }
00179 return 0;
00180 }
00181
00182 static int HighDegree(int node1[], int node2[], int nedges, int iwork[], int nnodes){
00183 int i,nmax=0,maxdegree=-1;
00184 for (i=0;i<nnodes;i++) iwork[i]=0;
00185 for (i=0;i<nedges;i++){
00186 iwork[node1[i]]++; iwork[node2[i]]++;
00187 }
00188 for (i=0;i<nnodes;i++){ if (iwork[i]>maxdegree){nmax=i; maxdegree=iwork[i];} }
00189 return nmax;
00190 }
00191
00192 static int First(int coloring[], int nnodes){
00193 int i,nmax=nnodes;
00194 for (i=0;i<nnodes;i++){
00195 if (coloring[i]==0){
00196 nmax=i; return nmax;
00197 }
00198 }
00199 return -1;
00200 }
00201
00202 static int RandomizedColor(DSDP dsdp, SDPCone sdpcone, int nnodes, int node1[], int node2[], int nedges){
00203 int i,j,nodek,nn,info,flag,coloring=0,maxvertex;
00204 int *degree,*color,*cgroup,ngsize,uncolored=nnodes;
00205 int tnedges=nedges;
00206 double *xrow,*xptr;
00207 orderVec *vorder;
00208
00209 xrow=(double*)malloc(nnodes*sizeof(double));
00210 color=(int*)malloc(nnodes*sizeof(int));
00211 cgroup=(int*)malloc(nnodes*sizeof(int));
00212 degree=(int*)malloc(nnodes*sizeof(int));
00213 vorder=(orderVec*)malloc(nnodes*sizeof(orderVec));
00214
00215 for (i=0;i<nnodes;i++){ color[i]=0;}
00216 info=DSDPComputeX(dsdp);
00217 info=SDPConeGetXArray(sdpcone,0,&xptr,&nn);
00218
00219 while (uncolored>0){
00220
00221 coloring++;
00222
00223 maxvertex=First(color,nnodes);
00224 maxvertex=HighDegree(node1,node2,tnedges,degree,nnodes);
00225
00226 cgroup[0]=maxvertex;ngsize=1;
00227
00228 info=GetXRow(xptr,xrow,maxvertex,nnodes);
00229
00230 for (i=0;i<nnodes;i++){vorder[i].index=i; vorder[i].val = xrow[i];}
00231 qsort( (void*)vorder, nnodes, sizeof(orderVec), cut_comp);
00232
00233 for (i=0;i<nnodes;i++){
00234 nodek=vorder[i].index;
00235 if (color[nodek]==0){
00236 for (flag=0,j=0;j<ngsize;j++){
00237 if (Connected(nodek,cgroup[j],node1,node2,tnedges) ){flag=1;}
00238 }
00239 if (flag==0){ cgroup[ngsize]=nodek; ngsize++; }
00240 }
00241 }
00242 for (i=0;i<ngsize;i++){
00243 color[cgroup[i]]=coloring; uncolored--;
00244 RemoveNode(cgroup[i],node1,node2,&tnedges);
00245 }
00246 }
00247 printf("\nCOLORS: %d\n",coloring);
00248 free(xrow);
00249 free(color);
00250 free(cgroup);
00251 free(degree);
00252 free(vorder);
00253 return(0);
00254 }
00255
00256
00257 #define BUFFERSIZ 100
00258
00259 #undef __FUNCT__
00260 #define __FUNCT__ "ParseGraphline"
00261 static int ParseGraphline(char * thisline,int *row,int *col,double *value,
00262 int *gotem){
00263
00264 int temp;
00265 int rtmp,coltmp;
00266
00267 rtmp=-1, coltmp=-1, *value=0.0;
00268 temp=sscanf(thisline,"%d %d %lf",&rtmp,&coltmp,value);
00269 if (temp==3 && rtmp>0 && coltmp>0) *gotem=3;
00270 else if (temp==2 && rtmp>0 && coltmp>0){ *value = 1.0; *gotem=3;}
00271 else *gotem=0;
00272 *row=rtmp-1; *col=coltmp-1;
00273
00274 return(0);
00275 }
00276
00277
00278 #undef __FUNCT__
00279 #define __FUNCT__ "ReadGraph"
00280 int ReadGraph(char* filename,int *nnodes, int *nedges,
00281 int**n1, int ** n2, double **wght){
00282
00283 FILE*fp;
00284 char thisline[BUFFERSIZ]="*";
00285 int i,k=0,line=0,nodes,edges,gotone=3;
00286 int *node1,*node2;
00287 double *weight;
00288 int info,row,col;
00289 double value;
00290
00291 fp=fopen(filename,"r");
00292 if (!fp){printf("Cannot open file %s !",filename);return(1);}
00293
00294 while(!feof(fp) && (thisline[0] == '*' || thisline[0] == '"') ){
00295 fgets(thisline,BUFFERSIZ,fp); line++;
00296 }
00297
00298 if (sscanf(thisline,"%d %d",&nodes, &edges)!=2){
00299 printf("First line must contain the number of nodes and number of edges\n");
00300 return 1;
00301 }
00302
00303 node1=(int*)malloc(edges*sizeof(int));
00304 node2=(int*)malloc(edges*sizeof(int));
00305 weight=(double*)malloc(edges*sizeof(double));
00306
00307 for (i=0; i<edges; i++){ node1[i]=0;node2[i]=0;weight[i]=0.0;}
00308
00309 while(!feof(fp) && gotone){
00310 thisline[0]='\0';
00311 fgets(thisline,BUFFERSIZ,fp); line++;
00312 info = ParseGraphline(thisline,&row,&col,&value,&gotone);
00313 if (gotone && value!=0.0 && k<edges &&
00314 col < nodes && row < nodes && col >= 0 && row >= 0){
00315 if (row<col){info=row;row=col;col=info;}
00316 if (row == col){}
00317 else {
00318 node1[k]=row; node2[k]=col;
00319 weight[k]=value; k++;
00320 }
00321 } else if (gotone && k>=edges) {
00322 printf("To many edges in file.\nLine %d, %s\n",line,thisline);
00323 return 1;
00324 } else if (gotone&&(col >= nodes || row >= nodes || col < 0 || row < 0)){
00325 printf("Invalid node number.\nLine %d, %s\n",line,thisline);
00326 return 1;
00327 }
00328 }
00329 *nnodes=nodes; *nedges=edges;
00330 *n1=node1; *n2=node2; *wght=weight;
00331 return 0;
00332 }