Logo Search packages:      
Sourcecode: dsdp version File versions  Download package

int MaxCut ( int  nnodes,
int  nedged,
int  node1[],
int  node2[],
double  weight[] 
)

Formulate and solve the SDP relaxation of the Maximum Cut problem.

Parameters:
nnodes number of nodes in graph
nedges number of edges in graph
node1 first node of each edge
node2 second node of each edge
weight weight of each edge
Note:
This routine is an example! It is not part of the solver library.

Definition at line 51 of file maxcut.c.

References DSDP_INFEASIBLE_START, DSDPComputeX(), DSDPCreate(), DSDPDestroy(), DSDPGetY(), DSDPReuseMatrix(), DSDPSetDualObjective(), DSDPSetGapTolerance(), DSDPSetPNormTolerance(), DSDPSetPotentialParameter(), DSDPSetR0(), DSDPSetStandardMonitor(), DSDPSetup(), DSDPSetY0(), DSDPSetZBar(), DSDPSolve(), DSDPStopReason(), MaxCutRandomized(), SDPConeAddASparseVecMat(), SDPConeGetXArray(), SDPConeSetASparseVecMat(), SDPConeSetBlockSize(), and SDPConeViewDataMatrix().

                                                                            {
  
  int i,info;
  int *indd,*iptr;
  double *yy,*val,*diag,tval=0;
  DSDPTerminationReason reason;
  SDPCone sdpcone;
  DSDP dsdp;
  
  info = DSDPCreate(nnodes,&dsdp); 
  info = DSDPCreateSDPCone(dsdp,1,&sdpcone);

  if (info){ printf("Out of memory\n"); return 1; }

  info = SDPConeSetBlockSize(sdpcone,0,nnodes);


  /* Formulate the problem from the data */
  /* 
     Diagonal elements equal 1.0 
     Create Constraint matrix A_i for i=1, ..., nnodes.
     that has a single nonzero element.
  */
  diag=(double*)malloc(nnodes*sizeof(double));
  iptr=(int*)malloc(nnodes*sizeof(int));
  for (i=0;i<nnodes;i++){
    iptr[i]=i*(i+1)/2+i; 
    diag[i]=1.0;
  }

  for (i=0;i<nnodes;i++){
    info = DSDPSetDualObjective(dsdp,i+1,1.0);
    info = SDPConeSetASparseVecMat(sdpcone,0,i+1,nnodes,1.0,0,iptr+i,diag+i,1);
    if (0==1){
      printf("Matrix: %d\n",i+1);
      info = SDPConeViewDataMatrix(sdpcone,0,i+1);
    }
  }

  /* C matrix is the Laplacian of the adjacency matrix */
  /* Also compute a feasible initial point y such that S>=0 */ 
  yy=(double*)malloc(nnodes*sizeof(double));
  for (i=0;i<nnodes;i++){yy[i]=0.0;}
  indd=(int*)malloc((nnodes+nedges)*sizeof(int));
  val=(double*)malloc((nnodes+nedges)*sizeof(double));
  for (i=0;i<nnodes+nedges;i++){indd[i]=0;}
  for (i=0;i<nnodes;i++){indd[nedges+i]=i*(i+1)/2+i;}
  for (i=0;i<nnodes+nedges;i++){val[i]=0;}
  for (i=0;i<nedges;i++){
    indd[i]=(node1[i])*(node1[i]+1)/2 + node2[i];
    tval+=fabs(weight[i]);
    val[i]=weight[i]/4;
    val[nedges+node1[i]]-=weight[i]/4;
    val[nedges+node2[i]]-=weight[i]/4;
    yy[node1[i]]-= fabs(weight[i]/2);
    yy[node2[i]]-= fabs(weight[i]/2);
  }

  if (0){
    info = SDPConeSetASparseVecMat(sdpcone,0,0,nnodes,1.0,0,indd,val,nedges+nnodes);
  } else { /* Equivalent */
    info = SDPConeSetASparseVecMat(sdpcone,0,0,nnodes,1.0,0,indd,val,nedges);
    info = SDPConeAddASparseVecMat(sdpcone,0,0,nnodes,1.0,0,indd+nedges,val+nedges,nnodes);
  } 
  if (0==1){ info = SDPConeViewDataMatrix(sdpcone,0,0);}
  
  /* Initial Point */
  info = DSDPSetR0(dsdp,0.0);
  info = DSDPSetZBar(dsdp,10*tval+1.0);
  for (i=0;i<nnodes; i++){ 
    info = DSDPSetY0(dsdp,i+1,10*yy[i]);
  }
  if (info) return info;

  /* Get read to go */
  info=DSDPSetGapTolerance(dsdp,0.001);
  info=DSDPSetPotentialParameter(dsdp,5);
  info=DSDPReuseMatrix(dsdp,0);
  info=DSDPSetPNormTolerance(dsdp,1.0);
  /*
  info = TCheckArgs(dsdp,argc,argv);
  */

  if (info){ printf("Out of memory\n"); return 1; }
  info = DSDPSetStandardMonitor(dsdp,1);

  info = DSDPSetup(dsdp);
  if (info){ printf("Out of memory\n"); return 1; }

  info = DSDPSolve(dsdp);
  if (info){ printf("Numerical error\n"); return 1; }
  info = DSDPStopReason(dsdp,&reason); 

  if (reason!=DSDP_INFEASIBLE_START){ /* Randomized solution strategy */
    info=MaxCutRandomized(sdpcone,nnodes);
    if (0==1){ /* Look at the solution */
      int n; double *xx,*y=diag;
      info=DSDPGetY(dsdp,y,nnodes);
      info=DSDPComputeX(dsdp);DSDPCHKERR(info);
      info=SDPConeGetXArray(sdpcone,0,&xx,&n);
    }
  }
  info = DSDPDestroy(dsdp);

  free(iptr);
  free(yy);
  free(indd);
  free(val);
  free(diag);
  
  return 0;
} /* main */


Generated by  Doxygen 1.6.0   Back to index