maxcut.c

Go to the documentation of this file.
00001 #include "dsdp5.h"
00002 #include <string.h>
00003 #include <math.h>
00004 #include <stdio.h>
00005 #include <stdlib.h>
00012 char help[]="\
00013 DSDP Usage: maxcut <graph filename> \n\
00014                    -gaptol <relative duality gap: default is 0.001> \n\
00015                   -maxit <maximum iterates: default is 200> \n";
00016 
00017 static int ReadGraph(char*,int *, int *,int**, int **, double **); 
00018 static int TCheckArgs(DSDP,int,char **);
00019 static int ParseGraphline(char*,int*,int*,double*,int*);
00020 int MaxCutRandomized(SDPCone,int);
00021 int MaxCut(int,int, int[], int[], double[]);
00022 
00023 
00024 int main(int argc,char *argv[]){
00025   int info;
00026   int *node1,*node2,nedges,nnodes;
00027   double *weight;
00028 
00029   if (argc<2){ printf("%s",help); return(1); }
00030 
00031   info = ReadGraph(argv[1],&nnodes,&nedges,&node1,&node2,&weight);
00032   if (info){ printf("Problem reading file\n"); return 1; }
00033 
00034   MaxCut(nnodes,nedges,node1,node2,weight);
00035 
00036   free(node1);free(node2);free(weight);
00037   return 0;
00038 }
00039 
00051 int MaxCut(int nnodes,int nedges, int node1[], int node2[], double weight[]){
00052   
00053   int i,info;
00054   int *indd,*iptr;
00055   double *yy,*val,*diag,tval=0;
00056   DSDPTerminationReason reason;
00057   SDPCone sdpcone;
00058   DSDP dsdp;
00059   
00060   info = DSDPCreate(nnodes,&dsdp); 
00061   info = DSDPCreateSDPCone(dsdp,1,&sdpcone);
00062 
00063   if (info){ printf("Out of memory\n"); return 1; }
00064 
00065   info = SDPConeSetBlockSize(sdpcone,0,nnodes);
00066 
00067 
00068   /* Formulate the problem from the data */
00069   /* 
00070      Diagonal elements equal 1.0 
00071      Create Constraint matrix A_i for i=1, ..., nnodes.
00072      that has a single nonzero element.
00073   */
00074   diag=(double*)malloc(nnodes*sizeof(double));
00075   iptr=(int*)malloc(nnodes*sizeof(int));
00076   for (i=0;i<nnodes;i++){
00077     iptr[i]=i*(i+1)/2+i; 
00078     diag[i]=1.0;
00079   }
00080 
00081   for (i=0;i<nnodes;i++){
00082     info = DSDPSetDualObjective(dsdp,i+1,1.0);
00083     info = SDPConeSetASparseVecMat(sdpcone,0,i+1,nnodes,1.0,0,iptr+i,diag+i,1);
00084     if (0==1){
00085       printf("Matrix: %d\n",i+1);
00086       info = SDPConeViewDataMatrix(sdpcone,0,i+1);
00087     }
00088   }
00089 
00090   /* C matrix is the Laplacian of the adjacency matrix */
00091   /* Also compute a feasible initial point y such that S>=0 */ 
00092   yy=(double*)malloc(nnodes*sizeof(double));
00093   for (i=0;i<nnodes;i++){yy[i]=0.0;}
00094   indd=(int*)malloc((nnodes+nedges)*sizeof(int));
00095   val=(double*)malloc((nnodes+nedges)*sizeof(double));
00096   for (i=0;i<nnodes+nedges;i++){indd[i]=0;}
00097   for (i=0;i<nnodes;i++){indd[nedges+i]=i*(i+1)/2+i;}
00098   for (i=0;i<nnodes+nedges;i++){val[i]=0;}
00099   for (i=0;i<nedges;i++){
00100     indd[i]=(node1[i])*(node1[i]+1)/2 + node2[i];
00101     tval+=fabs(weight[i]);
00102     val[i]=weight[i]/4;
00103     val[nedges+node1[i]]-=weight[i]/4;
00104     val[nedges+node2[i]]-=weight[i]/4;
00105     yy[node1[i]]-= fabs(weight[i]/2);
00106     yy[node2[i]]-= fabs(weight[i]/2);
00107   }
00108 
00109   if (0){
00110     info = SDPConeSetASparseVecMat(sdpcone,0,0,nnodes,1.0,0,indd,val,nedges+nnodes);
00111   } else { /* Equivalent */
00112     info = SDPConeSetASparseVecMat(sdpcone,0,0,nnodes,1.0,0,indd,val,nedges);
00113     info = SDPConeAddASparseVecMat(sdpcone,0,0,nnodes,1.0,0,indd+nedges,val+nedges,nnodes);
00114   } 
00115   if (0==1){ info = SDPConeViewDataMatrix(sdpcone,0,0);}
00116   
00117   /* Initial Point */
00118   info = DSDPSetR0(dsdp,0.0);
00119   info = DSDPSetZBar(dsdp,10*tval+1.0);
00120   for (i=0;i<nnodes; i++){ 
00121     info = DSDPSetY0(dsdp,i+1,10*yy[i]);
00122   }
00123   if (info) return info;
00124 
00125   /* Get read to go */
00126   info=DSDPSetGapTolerance(dsdp,0.001);
00127   info=DSDPSetPotentialParameter(dsdp,5);
00128   info=DSDPReuseMatrix(dsdp,0);
00129   info=DSDPSetPNormTolerance(dsdp,1.0);
00130   /*
00131   info = TCheckArgs(dsdp,argc,argv);
00132   */
00133 
00134   if (info){ printf("Out of memory\n"); return 1; }
00135   info = DSDPSetStandardMonitor(dsdp,1);
00136 
00137   info = DSDPSetup(dsdp);
00138   if (info){ printf("Out of memory\n"); return 1; }
00139 
00140   info = DSDPSolve(dsdp);
00141   if (info){ printf("Numerical error\n"); return 1; }
00142   info = DSDPStopReason(dsdp,&reason); 
00143 
00144   if (reason!=DSDP_INFEASIBLE_START){ /* Randomized solution strategy */
00145     info=MaxCutRandomized(sdpcone,nnodes);
00146     if (0==1){ /* Look at the solution */
00147       int n; double *xx,*y=diag;
00148       info=DSDPGetY(dsdp,y,nnodes);
00149       info=DSDPComputeX(dsdp);DSDPCHKERR(info);
00150       info=SDPConeGetXArray(sdpcone,0,&xx,&n);
00151     }
00152   }
00153   info = DSDPDestroy(dsdp);
00154 
00155   free(iptr);
00156   free(yy);
00157   free(indd);
00158   free(val);
00159   free(diag);
00160   
00161   return 0;
00162 } /* main */
00163 
00164 
00165 
00175 int MaxCutRandomized(SDPCone sdpcone,int nnodes){
00176   int i,j,derror,info;
00177   double dd,scal=RAND_MAX,*vv,*tt,*cc,ymin=0;
00178 
00179   vv=(double*)malloc(nnodes*sizeof(double));
00180   tt=(double*)malloc(nnodes*sizeof(double));
00181   cc=(double*)malloc((nnodes+2)*sizeof(double));
00182   info=SDPConeComputeXV(sdpcone,0,&derror);
00183   for (i=0;i<nnodes;i++){
00184       for (j=0;j<nnodes;j++){dd = (( rand())/scal - .5); vv[j] = tan(3.1415926*dd);}
00185       info=SDPConeXVMultiply(sdpcone,0,vv,tt,nnodes);
00186       for (j=0;j<nnodes;j++){if (tt[j]<0) tt[j]=-1; else tt[j]=1;}
00187       for (j=0;j<nnodes+2;j++){cc[j]=0;}
00188       info=SDPConeAddXVAV(sdpcone,0,tt,nnodes,cc,nnodes+2);
00189       if (cc[0]<ymin) ymin=cc[0];
00190   }
00191   printf("Best integer solution: %4.2f\n",ymin);
00192   free(vv); free(tt); free(cc);
00193 
00194   return(0);
00195 }
00196 
00197 static int TCheckArgs(DSDP dsdp,int nargs,char *runargs[]){
00198 
00199   int kk, info;
00200   
00201   info=DSDPSetOptions(dsdp,runargs,nargs);
00202   for (kk=1; kk<nargs-1; kk++){
00203     if (strncmp(runargs[kk],"-dloginfo",8)==0){
00204       info=DSDPLogInfoAllow(DSDP_TRUE,0);
00205     } else if (strncmp(runargs[kk],"-params",7)==0){
00206       info=DSDPReadOptions(dsdp,runargs[kk+1]);
00207     } else if (strncmp(runargs[kk],"-help",7)==0){
00208       printf("%s\n",help);
00209     } 
00210   }
00211 
00212   return 0;
00213 }
00214 
00215 
00216 #define BUFFERSIZ 100
00217 
00218 #undef __FUNCT__  
00219 #define __FUNCT__ "ParseGraphline"
00220 static int ParseGraphline(char * thisline,int *row,int *col,double *value, 
00221                           int *gotem){
00222 
00223   int temp;
00224   int rtmp,coltmp;
00225 
00226   rtmp=-1, coltmp=-1, *value=0.0;
00227   temp=sscanf(thisline,"%d %d %lf",&rtmp,&coltmp,value);
00228   if (temp==3 && rtmp>0 && coltmp>0) *gotem=3;
00229   else if (temp==2 && rtmp>0 && coltmp>0){ *value = 1.0; *gotem=3;}
00230   else *gotem=0;
00231   *row=rtmp-1; *col=coltmp-1;
00232 
00233   return(0);
00234 }
00235 
00236 
00237 #undef __FUNCT__  
00238 #define __FUNCT__ "ReadGraph"
00239 int ReadGraph(char* filename,int *nnodes, int *nedges,
00240               int**n1, int ** n2, double **wght){
00241 
00242   FILE*fp;
00243   char thisline[BUFFERSIZ]="*";
00244   int i,k=0,line=0,nodes,edges,gotone=3;
00245   int *node1,*node2;
00246   double *weight;
00247   int info,row,col;
00248   double value;
00249 
00250   fp=fopen(filename,"r");
00251   if (!fp){printf("Cannot open file %s !",filename);return(1);}
00252 
00253   while(!feof(fp) && (thisline[0] == '*' || thisline[0] == '"') ){
00254     fgets(thisline,BUFFERSIZ,fp); line++;
00255   }
00256 
00257   if (sscanf(thisline,"%d %d",&nodes, &edges)!=2){
00258     printf("First line must contain the number of nodes and number of edges\n");
00259     return 1;
00260   }
00261 
00262   node1=(int*)malloc(edges*sizeof(int));
00263   node2=(int*)malloc(edges*sizeof(int));
00264   weight=(double*)malloc(edges*sizeof(double));
00265 
00266   for (i=0; i<edges; i++){ node1[i]=0;node2[i]=0;weight[i]=0.0;}
00267   
00268   while(!feof(fp) && gotone){
00269     thisline[0]='\0';
00270     fgets(thisline,BUFFERSIZ,fp); line++;
00271     info = ParseGraphline(thisline,&row,&col,&value,&gotone);
00272     if (gotone && value!=0.0 && k<edges && 
00273         col < nodes && row < nodes && col >= 0 && row >= 0){
00274       if (row<col){info=row;row=col;col=info;}
00275       if (row == col){}
00276       else {
00277         node1[k]=row;        node2[k]=col;
00278         weight[k]=value;     k++;
00279       }
00280     } else if (gotone &&  k>=edges) {
00281       printf("To many edges in file.\nLine %d, %s\n",line,thisline);
00282       return 1;
00283     } else if (gotone&&(col >= nodes || row >= nodes || col < 0 || row < 0)){
00284       printf("Invalid node number.\nLine %d, %s\n",line,thisline);
00285       return 1;
00286     }
00287   }
00288   *nnodes=nodes; *nedges=edges;
00289   *n1=node1; *n2=node2; *wght=weight;
00290   return 0;
00291 }

Generated on Wed May 21 19:31:18 2008 for DSDP by  doxygen 1.5.5