12 Compute the Lovasz theta number for a graph. This number is an upper bound for \n\
13 the maximum clique of a graph, a lower bound for the mimimal graph coloring, and serves \n\
14 as a bound for several other combitorial graph problems. The number is the solution to \n\
15 a semidfinite program. \n\n\
16 The file should be the complement of the graph \n\
17 This file also demonstrates the use of customized data matrices in DSDP.\n\n\
18 DSDP Usage: theta <graph complement filename> \n";
24 #include "../src/sdp/dsdpdatamat_impl.h"
38 static int ReadGraphFromFile(
char*,
int*,
int*, EdgeMat*[]);
42 int main(
int argc,
char *argv[]){
58 int info,kk,nedges,nnodes;
63 if (argc<2){ printf(
"%s",help);
return(1); }
65 info = ReadGraphFromFile(argv[1],&nnodes,&nedges,&Edges);
66 if (info){ printf(
"Problem reading file\n");
return 1; }
69 info = DSDPCreateSDPCone(dsdp,1,&sdpcone);
73 if (info){ printf(
"Out of memory\n");
return 1; }
74 info =
SetThetaData(dsdp, sdpcone, nnodes, nedges, Edges);
75 if (info){ printf(
"Out of memory\n");
return 1; }
81 for (kk=1; kk<argc-1; kk++){
82 if (strncmp(argv[kk],
"-dloginfo",8)==0){
84 }
else if (strncmp(argv[kk],
"-params",7)==0){
86 }
else if (strncmp(argv[kk],
"-help",5)==0){
92 if (info){ printf(
"Out of memory\n");
return 1; }
94 if (info){ printf(
"Monitor Problem \n");
return 1; }
97 if (info){ printf(
"Out of memory\n");
return 1; }
101 if (info){ printf(
"Numberical error\n");
return 1; }
133 info=OneMatOpsInitialize(&onematops);
137 info=EdgeMatOperationsInitialize(&edgematop);
138 for (i=0; i<edges; i++){
142 if (info)
return info;
146 info = IdentitymatOperationsInitialize(&identitymatops);
149 info =
DSDPSetY0(dsdp,edges+1,-10*nodes-1);
157 #define BUFFERSIZ 100
160 #define __FUNCT__ "ParseGraphline"
161 static int ParseGraphline(
char * thisline,
int *row,
int *col,
double *value,
167 rtmp=-1, coltmp=-1, *value=0.0;
168 temp=sscanf(thisline,
"%d %d %lf",&rtmp,&coltmp,value);
169 if (temp==3 && rtmp>0 && coltmp>0) *gotem=3;
170 else if (temp==2 && rtmp>0 && coltmp>0){ *value = 1.0; *gotem=3;}
172 *row=rtmp-1; *col=coltmp-1;
173 if (*gotem && (*col < 0 || *row < 0)){
174 printf(
"Node Number must be positive.\n, %s\n",thisline);
180 #define __FUNCT__ "ReadGraphFromFile"
181 static int ReadGraphFromFile(
char* filename,
int *nnodes,
int *nedges, EdgeMat**EE){
184 char thisline[BUFFERSIZ]=
"*";
185 int i,k=0,line=0,nodes,edges,gotone=3;
190 fp=fopen(filename,
"r");
191 if (!fp){ printf(
"Cannot open file %s !",filename);
return(1); }
193 while(!feof(fp) && (thisline[0] ==
'*' || thisline[0] ==
'"') ){
194 fgets(thisline,BUFFERSIZ,fp); line++; }
196 if (sscanf(thisline,
"%d %d",&nodes, &edges)!=2){
197 printf(
"First line must contain the number of nodes and number of edges\n");
201 E=(EdgeMat*)malloc(edges*
sizeof(EdgeMat)); *EE=E;
202 for (i=0; i<edges; i++){ E[i].v1=0; E[i].v2=0; }
204 while(!feof(fp) && gotone){
206 fgets(thisline,BUFFERSIZ,fp); line++;
207 info = ParseGraphline(thisline,&row,&col,&value,&gotone);
208 if (gotone && k<edges &&
209 col < nodes && row < nodes && col >= 0 && row >= 0){
210 if (row > col){i=row;row=col;col=i;}
212 else { E[k].v1=row; E[k].v2=col; k++;}
213 }
else if (gotone && k>=edges) {
214 printf(
"To many edges in file.\nLine %d, %s\n",line,thisline);
216 }
else if (gotone&&(col >= nodes || row >= nodes || col < 0 || row < 0)){
217 printf(
"Invalid node number.\nLine %d, %s\n",line,thisline);
222 *nnodes=nodes; *nedges=k;
229 static int EdgeMatDestroy(
void*);
230 static int EdgeMatView(
void*);
231 static int EdgeMatVecVec(
void*,
double[],
int,
double *);
232 static int EdgeMatDot(
void*,
double[],
int,
int,
double *);
233 static int EdgeMatGetRank(
void*,
int*,
int);
234 static int EdgeMatFactor(
void*);
235 static int EdgeMatGetEig(
void*,
int,
double*,
double[],
int,
int[],
int*);
236 static int EdgeMatAddRowMultiple(
void*,
int,
double,
double[],
int);
237 static int EdgeMatAddMultiple(
void*,
double,
double[],
int,
int);
238 static int EdgeMatGetRowNnz(
void*,
int,
int[],
int*,
int);
240 static int EdgeMatDestroy(
void* AA){
243 static int EdgeMatVecVec(
void* A,
double x[],
int n,
double *v){
244 EdgeMat*E=(EdgeMat*)A;
245 *v=2*x[E->v1]*x[E->v2];
248 static int EdgeMatDot(
void* A,
double x[],
int nn,
int n,
double *v){
249 EdgeMat*E=(EdgeMat*)A;
250 int k=E->v2*(E->v2+1)/2 + E->v1;
254 static int EdgeMatView(
void* A){
255 EdgeMat*E=(EdgeMat*)A;
256 printf(
" Row: %d, Column: %d\n",E->v1,E->v2);
257 printf(
" Row: %d, Column: %d\n",E->v2,E->v1);
260 static int EdgeMatFactor(
void* A){
263 static int EdgeMatGetRank(
void *A,
int*rank,
int n){
267 static int EdgeMatGetEig(
void*A,
int neig,
double *eig,
double v[],
int n,
int spind[],
int *nind){
268 EdgeMat*E=(EdgeMat*)A;
269 double tt=1.0/sqrt(2.0);
270 memset((
void*)v,0,(n)*
sizeof(
double));
271 memset((
void*)spind,0,(n)*
sizeof(
int));
273 v[E->v1]=tt;v[E->v2]=tt;*eig=1;
274 spind[0]=E->v1;spind[1]=E->v2; *nind=2;
276 v[E->v1]=-tt;v[E->v2]=tt;*eig=-1;
277 spind[0]=E->v1;spind[1]=E->v2; *nind=2;
278 }
else { *eig=0;*nind=0;}
281 static int EdgeMatGetRowNnz(
void*A,
int nrow,
int nz[],
int *nnzz,
int n){
282 EdgeMat*E=(EdgeMat*)A;
283 if (nrow==E->v1){ nz[E->v2]++; *nnzz=1;}
284 else if (nrow==E->v2){nz[E->v1]++; *nnzz=1;}
288 static int EdgeMatAddRowMultiple(
void*A,
int nrow,
double dd,
double rrow[],
int n){
289 EdgeMat*E=(EdgeMat*)A;
290 if (nrow==E->v1){ rrow[E->v2]+=dd;}
291 else if (nrow==E->v2){rrow[E->v1]+=dd;}
294 static int EdgeMatAddMultiple(
void*A,
double dd,
double vv[],
int nn,
int n){
295 EdgeMat*E=(EdgeMat*)A;
296 int k=E->v2*(E->v2+1)/2 + E->v1;
300 static int EdgeMatFNorm(
void*A,
int n,
double *fnorm){
304 static int EdgeMatCountNonzeros(
void*A,
int *nnz,
int n){
308 static const char *datamatname=
"THETA EDGE MATRIX";
309 static int EdgeMatOperationsInitialize(
struct DSDPDataMat_Ops* edgematoperator){
311 if (edgematoperator==NULL)
return 0;
313 edgematoperator->matfactor1=EdgeMatFactor;
314 edgematoperator->matgetrank=EdgeMatGetRank;
315 edgematoperator->matgeteig=EdgeMatGetEig;
316 edgematoperator->matvecvec=EdgeMatVecVec;
317 edgematoperator->matrownz=EdgeMatGetRowNnz;
318 edgematoperator->matdot=EdgeMatDot;
319 edgematoperator->matfnorm2=EdgeMatFNorm;
320 edgematoperator->matnnz=EdgeMatCountNonzeros;
321 edgematoperator->mataddrowmultiple=EdgeMatAddRowMultiple;
322 edgematoperator->mataddallmultiple=EdgeMatAddMultiple;
323 edgematoperator->matdestroy=EdgeMatDestroy;
324 edgematoperator->matview=EdgeMatView;
325 edgematoperator->matname=datamatname;
326 edgematoperator->id=25;
331 static int OneMatDestroy(
void*);
332 static int OneMatView(
void*);
333 static int OneMatVecVec(
void*,
double[],
int,
double *);
334 static int OneMatDot(
void*,
double[],
int,
int,
double *);
335 static int OneMatGetRank(
void*,
int*,
int);
336 static int OneMatFactor(
void*);
337 static int OneMatGetEig(
void*,
int,
double*,
double[],
int,
int[],
int*);
338 static int OneMatRowNnz(
void*,
int,
int[],
int*,
int);
339 static int OneMatAddRowMultiple(
void*,
int,
double,
double[],
int);
340 static int OneMatAddMultiple(
void*,
double,
double[],
int,
int);
343 static int OneMatFactor(
void*A){
return 0;}
344 static int OneMatGetRank(
void *A,
int *rank,
int n){*rank=1;
return 0;}
345 static int OneMatFNorm2(
void*AA,
int n,
double *fnorm2){*fnorm2=1.0*n*n;
return 0;}
346 static int OneMatCountNonzeros(
void*AA,
int *nnz,
int n){*nnz=n*n;
return 0;}
347 static int OneMatDot(
void* A,
double x[],
int nn,
int n,
double *v){
351 for (j=0;j<i;j++,x++){dtmp+= (*x);}
357 static int OneMatVecVec(
void* A,
double x[],
int n,
double *v){
360 for (i=0; i<n; i++){dtmp+=x[i];}
364 static int OneMatAddMultiple(
void*A,
double ddd,
double vv[],
int nn,
int n){
367 for (j=0;j<i;j++,vv++){(*vv)+=-ddd;}
372 static int OneMatAddRowMultiple(
void*A,
int nrow,
double ddd,
double row[],
int n){
374 for (i=0;i<n;i++){row[i] -= ddd;}
378 static int OneMatGetEig(
void*A,
int neig,
double *eig,
double v[],
int n,
int spind[],
int *nind){
380 if (neig==0){ *eig=-1;
for (i=0;i<n;i++){ v[i]=1.0; spind[i]=i;} *nind=n;
381 }
else { *eig=0;
for (i=0;i<n;i++){ v[i]=0.0; } *nind=0;
385 static int OneMatRowNnz(
void*A,
int row,
int nz[],
int *nnz,
int n){
387 for (i=0;i<n;i++){ nz[i]++; }
391 static int OneMatView(
void* AA){
392 printf(
"Every element of the matrix is the same: -1\n");
395 static int OneMatDestroy(
void* A){
return 0;}
397 static const char *mat1name=
"THETA ALL ELEMENTS EQUAL -ONE";
400 if (mat1ops==NULL)
return 0;
402 mat1ops->matfactor1=OneMatFactor;
403 mat1ops->matgetrank=OneMatGetRank;
404 mat1ops->matgeteig=OneMatGetEig;
405 mat1ops->matvecvec=OneMatVecVec;
406 mat1ops->matdot=OneMatDot;
407 mat1ops->mataddrowmultiple=OneMatAddRowMultiple;
408 mat1ops->mataddallmultiple=OneMatAddMultiple;
409 mat1ops->matdestroy=OneMatDestroy;
410 mat1ops->matview=OneMatView;
411 mat1ops->matrownz=OneMatRowNnz;
412 mat1ops->matfnorm2=OneMatFNorm2;
413 mat1ops->matnnz=OneMatCountNonzeros;
415 mat1ops->matname=mat1name;
420 static int IdentityMatDestroy(
void*);
421 static int IdentityMatView(
void*);
422 static int IdentityMatVecVec(
void*,
double[],
int,
double *);
423 static int IdentityMatDot(
void*,
double[],
int,
int,
double *);
424 static int IdentityMatGetRank(
void*,
int*,
int);
425 static int IdentityMatFactor(
void*);
426 static int IdentityMatGetEig(
void*,
int,
double*,
double[],
int,
int[],
int*);
427 static int IdentityMatAddRowMultiple(
void*,
int,
double,
double[],
int);
428 static int IdentityMatAddMultiple(
void*,
double,
double[],
int,
int);
429 static int IdentityMatGetRowNnz(
void*,
int,
int[],
int*,
int);
431 static int IdentityMatDestroy(
void* AA){
return 0;}
432 static int IdentityMatFNorm2(
void* AA,
int n,
double *v){*v=1.0*n;
return 0;}
433 static int IdentityMatGetRank(
void *AA,
int*rank,
int n){*rank=n;
return 0;}
434 static int IdentityMatFactor(
void*A){
return 0;}
435 static int IdentityMatCountNonzeros(
void*A,
int *nnz,
int n){*nnz=n;
return 0;}
436 static int IdentityMatVecVec(
void* AA,
double x[],
int n,
double *v){
439 for (i=0;i<n;i++){ *v+=x[i]*x[i]; }
442 static int IdentityMatDot(
void* AA,
double x[],
int nn,
int n,
double *v){
445 for (i=0;i<n;i++){ vv+=x[((i+1)*(i+2))/2-1];}
449 static int IdentityMatView(
void* AA){
450 printf(
"Identity matrix: All Diagonal elements equal 1.0\n");
453 static int IdentityMatGetEig(
void*AA,
int neig,
double *eig,
double v[],
int n,
int spind[],
int *nind){
454 if (neig<0 || neig>=n){ *eig=0;
return 0;}
455 memset((
void*)v,0,(n)*
sizeof(
double));
456 *eig=1.0; v[neig]=1.0; spind[0]=neig; *nind=1;
459 static int IdentityMatGetRowNnz(
void*A,
int nrow,
int nz[],
int *nnzz,
int n){
460 if (nrow>=0 && nrow < n){
466 static int IdentityMatAddRowMultiple(
void*A,
int nrow,
double dd,
double rrow[],
int n){
467 rrow[nrow] += dd;
return 0;
469 static int IdentityMatAddMultiple(
void*A,
double dd,
double vv[],
int nn,
int n){
471 for (i=0;i<n;i++){ vv[(i+1)*(i+2)/2-1] += dd;}
475 static const char *eyematname=
"THETA IDENTITY MATRIX";
476 static int IdentitymatOperationsInitialize(
struct DSDPDataMat_Ops* imatops){
478 if (imatops==NULL)
return 0;
480 imatops->matfactor1=IdentityMatFactor;
481 imatops->matgetrank=IdentityMatGetRank;
482 imatops->matgeteig=IdentityMatGetEig;
483 imatops->matvecvec=IdentityMatVecVec;
484 imatops->matrownz=IdentityMatGetRowNnz;
485 imatops->matdot=IdentityMatDot;
486 imatops->matfnorm2=IdentityMatFNorm2;
487 imatops->matnnz=IdentityMatCountNonzeros;
488 imatops->mataddrowmultiple=IdentityMatAddRowMultiple;
489 imatops->mataddallmultiple=IdentityMatAddMultiple;
490 imatops->matdestroy=IdentityMatDestroy;
491 imatops->matview=IdentityMatView;
493 imatops->matname=eyematname;