DSDP
|
00001 #include "dsdp5.h" 00002 #include <stdio.h> 00003 #include <string.h> 00004 #include <math.h> 00005 #include <stdlib.h> 00006 00011 char help[]="\n\ 00012 Compute the stable set for a graph. \n\n\ 00013 DSDP Usage: stable <graph filename> \n"; 00014 00015 typedef struct { 00016 double v[3]; 00017 int indd[3]; 00018 } EdgeMat; 00019 00020 static int ReadGraphFromFile(char*,int*, int*, EdgeMat*[]); 00021 int SetStableSetData(DSDP, SDPCone, int, int, EdgeMat[]); 00022 int StableSet(int argc,char *argv[]); 00023 int StableRandomized(SDPCone,int, int, EdgeMat[]); 00024 #define CHKERR(a) { if (a){printf("DSDP Numerical Error or Memory Problem"); return 2;} } 00025 00026 int main(int argc,char *argv[]){ 00027 int info; 00028 info=StableSet(argc,argv); 00029 return 0; 00030 } 00031 00040 int StableSet(int argc,char *argv[]){ 00041 00042 int info,kk,nedges,nnodes; 00043 EdgeMat *Edges; 00044 DSDP dsdp; 00045 SDPCone sdpcone; 00046 00047 if (argc<2){ printf("%s",help); return(1); } 00048 00049 info = ReadGraphFromFile(argv[1],&nnodes,&nedges,&Edges); 00050 if (info){ printf("Problem reading file\n"); return 1; } 00051 00052 info = DSDPCreate(nedges+nnodes+1,&dsdp);CHKERR(info); 00053 info = DSDPCreateSDPCone(dsdp,1,&sdpcone);CHKERR(info); 00054 info = SDPConeSetBlockSize(sdpcone,0,nnodes+1);CHKERR(info); 00055 info = SDPConeUsePackedFormat(sdpcone,0);CHKERR(info); 00056 info = SDPConeSetSparsity(sdpcone,0,nedges+nnodes+1);CHKERR(info); 00057 info = SetStableSetData(dsdp, sdpcone, nnodes, nedges, Edges); 00058 if (info){ printf("Problem setting data\n"); return 1; } 00059 00060 info = DSDPSetGapTolerance(dsdp,0.0001);CHKERR(info); 00061 info = DSDPSetZBar(dsdp,1e10*nnodes+1);CHKERR(info); 00062 info = DSDPReuseMatrix(dsdp,10);CHKERR(info); 00063 00064 for (kk=1; kk<argc-1; kk++){ 00065 if (strncmp(argv[kk],"-dloginfo",8)==0){ 00066 info=DSDPLogInfoAllow(DSDP_TRUE,0);CHKERR(info); 00067 } else if (strncmp(argv[kk],"-params",7)==0){ 00068 info=DSDPReadOptions(dsdp,argv[kk+1]);CHKERR(info); 00069 } else if (strncmp(argv[kk],"-help",5)==0){ 00070 printf("%s\n",help); 00071 } 00072 } 00073 info=DSDPSetOptions(dsdp,argv,argc);CHKERR(info); 00074 00075 info = DSDPSetStandardMonitor(dsdp,1);CHKERR(info); 00076 00077 info = DSDPSetup(dsdp);CHKERR(info); 00078 if (0==1){info=SDPConeCheckData(sdpcone);} 00079 00080 info=DSDPSolve(dsdp); CHKERR(info); 00081 info=StableRandomized(sdpcone,nnodes,nedges,Edges); 00082 00083 info=DSDPComputeX(dsdp);CHKERR(info); 00084 00085 if (0==1){ /* Look at the solution */ 00086 int n; double *xx; 00087 info=SDPConeGetXArray(sdpcone,0,&xx,&n);CHKERR(info); 00088 } 00089 00090 info = DSDPDestroy(dsdp);CHKERR(info); 00091 free(Edges); 00092 00093 return 0; 00094 } /* main */ 00095 00107 int SetStableSetData(DSDP dsdp, SDPCone sdpcone, int nodes, int edges, EdgeMat Edge[]){ 00108 00109 int i,ii,info,nnodes=nodes+1; 00110 int *iptr,*iptr2; 00111 double *diag; 00112 00113 /* The c matrix has all elements equal to 1.0 */ 00114 diag=(double*)malloc(nnodes*sizeof(double)); 00115 iptr=(int*)malloc(nnodes*sizeof(int)); 00116 iptr2=(int*)malloc(nnodes*sizeof(int)); 00117 00118 ii=nodes*(nodes+1)/2; 00119 for (i=0;i<nnodes;i++){ 00120 diag[i]=1.0; 00121 iptr[i]=i*(i+1)/2+i; 00122 iptr2[i]=i; 00123 } 00124 info = SDPConeSetASparseVecMat(sdpcone,0,0,nnodes,-0.50,0,iptr,diag,nodes);CHKERR(info); 00125 info = SDPConeAddASparseVecMat(sdpcone,0,0,nnodes,-0.25,-ii,iptr2,diag,nodes);CHKERR(info); 00126 if (0){info=SDPConeViewDataMatrix(sdpcone,0,0);} 00127 /* Diagonal elements must equal 1 */ 00128 for (i=0;i<nnodes;i++){ 00129 info = DSDPSetDualObjective(dsdp,i+1,1.0);CHKERR(info); 00130 info = DSDPSetY0(dsdp,i+1,0.0);CHKERR(info); 00131 info = SDPConeSetASparseVecMat(sdpcone,0,i+1,nnodes,1.0,0,iptr+i,diag+i,1);CHKERR(info); 00132 } 00133 00134 /* 00135 For each edge connecting nodes i and j, 00136 X(i,i)+X(j,j)+X(i,j)+X(j,i)+X(i,n)+X(n,i)+X(j,n)+X(n,j)+X(n,n) = 1 00137 where nodes i,j numbered 0 ... n-1. 00138 */ 00139 for (i=0; i<edges; i++){ 00140 info = SDPConeAddARankOneMat(sdpcone,0,i+nnodes+1,nnodes,1.0,0,Edge[i].indd,Edge[i].v,3); 00141 if (0==1){info = SDPConeViewDataMatrix(sdpcone,0,i+nnodes+1);CHKERR(info);} 00142 info = DSDPSetDualObjective(dsdp,i+nnodes+1,1.0);CHKERR(info); 00143 info = DSDPSetY0(dsdp,i+nnodes+1,0.0);CHKERR(info); 00144 } 00145 00146 /* The initial point y is feasible and near the central path */ 00147 /* 00148 info = DSDPSetR0(dsdp,0); 00149 */ 00150 return(0); 00151 } 00152 00164 int StableRandomized(SDPCone sdpcone,int nodes, int edges, EdgeMat Edge[]){ 00165 int i,j,derror,info,nnodes=nodes+1; 00166 double dd,scal=RAND_MAX,*vv,*tt,*cc,ymin=0; 00167 int e0,e1,e2; 00168 00169 vv=(double*)malloc(nnodes*sizeof(double)); 00170 tt=(double*)malloc(nnodes*sizeof(double)); 00171 cc=(double*)malloc((edges+nnodes+2)*sizeof(double)); 00172 info=SDPConeComputeXV(sdpcone,0,&derror); 00173 for (i=0;i<nnodes;i++){ 00174 for (j=0;j<nnodes;j++){dd = (( rand())/scal - .5); vv[j] = tan(3.1415926*dd);} 00175 info=SDPConeXVMultiply(sdpcone,0,vv,tt,nnodes); 00176 for (j=0; j<edges; j++){ 00177 e0=Edge[j].indd[0];e1=Edge[j].indd[1];e2=Edge[j].indd[2]; 00178 if (tt[e0] * tt[e2] > 0 && tt[e1]*tt[e2] >0){ 00179 if ( fabs(tt[e0]-tt[e2]) > fabs(tt[e1]-tt[e2]) ){ 00180 tt[e0]*=-1; 00181 } else { 00182 tt[e1]*=-1; 00183 } 00184 } 00185 } 00186 for (j=0;j<nnodes;j++){if (tt[j]<0) tt[j]=-1; else tt[j]=1;} 00187 for (j=0;j<edges+nodes+1;j++){cc[j]=0;} 00188 info=SDPConeAddXVAV(sdpcone,0,tt,nnodes,cc,edges+nnodes+2); 00189 if (cc[0]<ymin) ymin=cc[0]; 00190 } 00191 printf("Stable Set Size: %4.0f\n",-ymin); 00192 free(vv); free(tt); free(cc); 00193 00194 return(0); 00195 } 00196 00197 00198 #define BUFFERSIZ 100 00199 00200 #undef __FUNCT__ 00201 #define __FUNCT__ "ParseGraphline" 00202 static int ParseGraphline(char * thisline,int *row,int *col,double *value, 00203 int *gotem){ 00204 00205 int temp; 00206 int rtmp,coltmp; 00207 00208 rtmp=-1, coltmp=-1, *value=0.0; 00209 temp=sscanf(thisline,"%d %d %lf",&rtmp,&coltmp,value); 00210 if (temp==3 && rtmp>0 && coltmp>0) *gotem=3; 00211 else if (temp==2 && rtmp>0 && coltmp>0){ *value = 1.0; *gotem=3;} 00212 else *gotem=0; 00213 *row=rtmp-1; *col=coltmp-1; 00214 if (*gotem && (*col < 0 || *row < 0)){ 00215 printf("Node Number must be positive.\n, %s\n",thisline); 00216 } 00217 return(0); 00218 } 00219 00220 #undef __FUNCT__ 00221 #define __FUNCT__ "ReadGraphFromFile" 00222 static int ReadGraphFromFile(char* filename,int *nnodes, int *nedges, EdgeMat**EE){ 00223 00224 FILE*fp; 00225 char thisline[BUFFERSIZ]="*"; 00226 int i,k=0,line=0,nodes,edges,gotone=3; 00227 int info,row,col; 00228 double value; 00229 EdgeMat *E; 00230 00231 fp=fopen(filename,"r"); 00232 if (!fp){ printf("Cannot open file %s !",filename); return(1); } 00233 00234 while(!feof(fp) && (thisline[0] == '*' || thisline[0] == '"') ){ 00235 fgets(thisline,BUFFERSIZ,fp); line++; } 00236 00237 if (sscanf(thisline,"%d %d",&nodes, &edges)!=2){ 00238 printf("First line must contain the number of nodes and number of edges\n"); 00239 return 1; 00240 } 00241 00242 E=(EdgeMat*)malloc(edges*sizeof(EdgeMat)); *EE=E; 00243 for (i=0; i<edges; i++){ 00244 E[i].v[0]=0; E[i].v[1]=0; E[i].v[2]=0; 00245 E[i].indd[0]=0; E[i].indd[1]=0; E[i].indd[2]=0; 00246 } 00247 00248 while(!feof(fp) && gotone){ 00249 thisline[0]='\0'; 00250 fgets(thisline,BUFFERSIZ,fp); line++; 00251 info = ParseGraphline(thisline,&row,&col,&value,&gotone); 00252 if (gotone && k<edges && 00253 col < nodes && row < nodes && col >= 0 && row >= 0){ 00254 if (row > col){i=row;row=col;col=i;} 00255 if (row == col){} 00256 else { 00257 E[k].indd[0]=row; E[k].indd[1]=col; E[k].indd[2]=nodes; 00258 E[k].v[0]=1.0; E[k].v[1]=1.0; E[k].v[2]=1.0; 00259 k++; 00260 } 00261 } else if (gotone && k>=edges) { 00262 printf("To many edges in file.\nLine %d, %s\n",line,thisline); 00263 return 1; 00264 } else if (gotone&&(col >= nodes || row >= nodes || col < 0 || row < 0)){ 00265 printf("Invalid node number.\nLine %d, %s\n",line,thisline); 00266 return 1; 00267 } 00268 } 00269 00270 *nnodes=nodes; *nedges=k; 00271 00272 return 0; 00273 } 00274