DSDP
dsdp.c
Go to the documentation of this file.
1 #include "mex.h"
2 #include "dsdp5.h"
3 #include <math.h>
4 
5 /* #define K_IN prhs[0] */
6 #define CA_IN prhs[1]
7 #define B_IN prhs[0]
8 #define PARS_IN prhs[2]
9 #define Y_IN prhs[3]
10 #define STAT_OUT plhs[0]
11 #define Y_OUT plhs[1]
12 #define X_OUT plhs[2]
13 static int DSDPPrintStats2(DSDP, void*);
14 static int CountNonzeroMatrices(double*,int,int,int,int*,int*,int*);
15 static int CheckForConstantMat(double*,int,int);
16 
17 static int printlevel=10;
18 #define CHKERR(a) { if (a){mexWarnMsgTxt("DSDP Numerical Error"); } }
19 
20 #define mexErrMsgTxt(a); mexPrintf("Error: "); mexWarnMsgTxt(a); return;
21 
35 void mexFunction(int nlhs, mxArray *plhs[],
36  int nrhs, const mxArray *prhs[]){
37 
38  mxArray *CA_cell_pr,*X_cell_pr;
39  const mxArray *OPTIONS_FIELD;
40  mxArray *STAT_FIELD;
41  int i,j,k,ii,itmp,index,info;
42  int *air, *ajc, *air2, *ajc2, *str,*iptr;
43  int nvars,nb,mC,nC,m1,n1,m2,n2;
44  int nsubs=2, subs[2];
45  int iscellCA;
46  int its,reuse=4,print_info=1,printsummary=0;
47  int ijnnz,spot,n,nn=0,nzmats,vecn;
48  int it1,it2;
49  int nsdpblocks=0,sdpblockj=0,sdpnmax=1,lpnmax=1,stat1=1,xmaker=0;
50  int sspot,nsubblocks,blockj;
51  int jj,tnnz,tnnz2;
52  int maxit=1000,fastblas=1,rpos=0,drho=1,iloginfo=0,aggressive=0;
53  double penalty=1e8,rho=4,zbar=1e10,cc=0,r0=-1,mu0=-1,ylow,yhigh,gaptol=1e-6,pnormtol=1e30;
54  double maxtrust=1e30,steptol=0.01,inftol=1e-8,lpb=1.0,dbound=1e20,infptol=1e-4;
55  double dtmp,pstep,dstep,pnorm,mu;
56  double *blockn,datanorm[3];
57  double *aval,*aval2,*bval,*yout,*y0=0,*xout,*stat;
58  double pobj,dobj,dinf;
59  DSDP dsdp;
60  SDPCone sdpcone=0;
61  DSDPTerminationReason reason;
62  DSDPSolutionType pdfeasible;
63  LPCone lpcone=0;
64  BCone bcone=0;
65  char conetype[30];
66  int nfields=25;
67  const char *fnames[25]={"stype","obj","pobj","dobj","stopcode","termcode","iterates","r","mu",
68  "pstep","dstep","pnorm","gaphist","infeashist","errors",
69  "datanorm","ynorm","boundy","penalty","tracex","reuse","rho","xy","xdy","xmu"};
70 
71  if (nrhs < 2){
72  mexErrMsgTxt("Two input arguments required. See help for details. ");}
73  if (nrhs > 4){
74  mexErrMsgTxt("Fewer input arguments required. See help for details. ");}
75  if (nlhs < 2){
76  mexErrMsgTxt("Two output arguments required. See help for details. ");}
77  if (nlhs > 3){
78  mexErrMsgTxt("Fewer output arguments required. See help for details. ");}
79 
80  if (!mxIsDouble(B_IN) || mxIsSparse(B_IN)){
81  mexErrMsgTxt("DSDP: 1ST input must be a dense vector of doubles"); }
82  nvars = mxGetM(B_IN);
83  nb = mxGetN(B_IN);
84  if (nb > 1){
85  mexErrMsgTxt("DESP: 1ST input must be a column vector"); }
86 
87  iscellCA = mxIsCell(CA_IN);
88  if (!iscellCA){
89  mexErrMsgTxt("DSDP: 2ND input must be a cell array"); }
90  mC = mxGetM(CA_IN);
91  nC = mxGetN(CA_IN);
92  if (nC != 3){
93  mexErrMsgTxt("DSDP: dimension of 2ND cell array is p x 3");}
94 
95  if (nrhs >2){
96  if(!mxIsStruct(PARS_IN)){
97  mexErrMsgTxt("3RD input `OPTIONS' should be a structure.");
98  }
99  }
100 
101  if (nrhs>3){
102  if (!mxIsDouble(Y_IN) || mxIsSparse(Y_IN)){
103  mexErrMsgTxt("DSDP: 4TH input must be a dense vector of doubles"); }
104  m1 = mxGetM(Y_IN);
105  n1 = mxGetN(Y_IN);
106  if (m1 != nvars || n1 != nb){
107  mexErrMsgTxt("DSDP: dimensions of 1ST and 4TH input not compatible");}
108  y0=mxGetPr(Y_IN);
109  if (!y0){
110  mexErrMsgTxt("DSDP: Cannot read 4TH argument");}
111  } else y0=0;
112 
113 
114  /* Check data */
115  for (j=0; j<mC; j++){
116  subs[0] = j; subs[1] = 0;
117  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
118  CA_cell_pr = mxGetCell(CA_IN,index);
119  if (!CA_cell_pr){
120  mexPrintf("??? DSDP DATA ERROR: Cell: %d, Column: %d\n",j+1,1);
121  mexErrMsgTxt("DSDP: Empty Cell. Missing String"); }
122  if (!mxIsChar(CA_cell_pr)){
123  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
124  mexErrMsgTxt("DSDP: First column of cells in 2ND argument are a string to determine cone type."); }
125  mxGetString(CA_cell_pr,conetype,20);
126 
127  subs[0] = j; subs[1] = 1;
128  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
129  CA_cell_pr = mxGetCell(CA_IN,index);
130  if (!CA_cell_pr){
131  mexPrintf("??? DSDP DATA ERROR: Cell: %d, Column: %d\n",j+1,2);
132  mexErrMsgTxt("DSDP: Empty Cell. Provide dimension of block."); }
133  if (mxIsSparse(CA_cell_pr)){
134  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
135  mexErrMsgTxt("DSDP: Second column in 2ND argument must be a dense array of scalars that specify dimension."); }
136  if (!mxIsDouble(CA_cell_pr)){
137  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
138  mexErrMsgTxt("DSDP: Second column in 2ND argument must specify dimension."); }
139  aval=mxGetPr(CA_cell_pr);
140 
141  subs[0] = j; subs[1] = 2;
142  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
143  CA_cell_pr = mxGetCell(CA_IN,index);
144  if (!CA_cell_pr){
145  mexPrintf("??? DSDP DATA ERROR: Cell: %d, Column: %d\n",j+1,3);
146  mexErrMsgTxt("DSDP: Empty Cell. Provide sparse data matrix."); }
147  if (!mxIsSparse(CA_cell_pr) || !mxIsDouble(CA_cell_pr)){
148  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
149  mexErrMsgTxt("DSDP: Third column in 2ND argument must be a real sparse data matrix."); }
150 
151  if (strcmp(conetype,"SDP")==0){
152  subs[0] = j; subs[1] = 1;
153  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
154  CA_cell_pr = mxGetCell(CA_IN,index);
155  aval=mxGetPr(CA_cell_pr);
156  it1 = mxGetM(CA_cell_pr);
157  it2 = mxGetN(CA_cell_pr);
158  if (it1!=1 && it2!=1){
159  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
160  mexErrMsgTxt("DSDP: Use a dense row vector in the second column in 2ND of the argument.");}
161 
162  subs[0] = j; subs[1] = 2;
163  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
164  CA_cell_pr = mxGetCell(CA_IN,index);
165  m1 = mxGetN(CA_cell_pr)-1;
166  if (m1 != nvars){
167  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
168  mexErrMsgTxt("DSDP: The matrix in the third column in 2ND argument must have number of columns equal to number of variables+1.");}
169  vecn = mxGetM(CA_cell_pr);
170  for (tnnz=0,i=0;i<it1*it2;i++){n=(int)aval[i]; tnnz=tnnz+n*(n+1)/2;}
171  if ( tnnz != vecn){
172  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
173  mexErrMsgTxt("DSDP: Check Dimensions: The columns of A and C cannot be converted into square matrices");}
174  nsdpblocks=nsdpblocks+it1*it2;
175  } else if (strcmp(conetype,"LP")==0){
176 
177  subs[0] = j; subs[1] = 1;
178  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
179  CA_cell_pr = mxGetCell(CA_IN,index);
180  it1 = mxGetM(CA_cell_pr);
181  it2 = mxGetN(CA_cell_pr);
182  /* THe sum of the dimensions should equal the number of constraints */
183 
184  subs[0] = j; subs[1] = 2;
185  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
186  CA_cell_pr = mxGetCell(CA_IN,index);
187  if (!mxIsSparse(CA_cell_pr)){
188  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
189  mexErrMsgTxt("DSDP: Matrices in the third column in 2ND argument must be sparse."); }
190  m1 = mxGetN(CA_cell_pr)-1;
191  if (m1 != nvars){
192  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
193  mexErrMsgTxt("DSDP: The matrix in the third column in 2ND argument must have number of columns equal to number of variables+1.");}
194  } else if (strcmp(conetype,"LB")==0 || strcmp(conetype,"UB")==0 ){
195 
196  subs[0] = j; subs[1] = 2;
197  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
198  CA_cell_pr = mxGetCell(CA_IN,index);
199  if (mxIsSparse(CA_cell_pr)){
200  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
201  mexErrMsgTxt("DSDP: Row vector in the third column in 2ND argument must be full."); }
202  it1 = mxGetM(CA_cell_pr);
203  it2 = mxGetN(CA_cell_pr);
204  if (it2 != 1){
205  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
206  mexErrMsgTxt("DSDP: The matrix in the third column in 2ND argument must have a single column of bounds.");}
207  if (it1 != nvars){
208  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
209  mexErrMsgTxt("DSDP: The column matrix in the third column in 2ND argument must contain a bound for each y variable.");}
210 
211  } else if (strcmp(conetype,"FIXED")==0){
212  int dim1;
213  subs[0] = j; subs[1] = 1;
214  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
215  CA_cell_pr = mxGetCell(CA_IN,index);
216  if (mxIsSparse(CA_cell_pr)){
217  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
218  mexErrMsgTxt("DSDP: Vector in the third column in 2ND argument must be full."); }
219  it1 = mxGetM(CA_cell_pr);
220  it2 = mxGetN(CA_cell_pr);
221  dim1=it1*it2;
222  if (it1 != 1){
223  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
224  mexErrMsgTxt("DSDP: Third column in 2ND argument must have 1 row,");}
225  subs[0] = j; subs[1] = 2;
226  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
227  CA_cell_pr = mxGetCell(CA_IN,index);
228  if (mxIsSparse(CA_cell_pr)){
229  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
230  mexErrMsgTxt("DSDP: Vector in the third column in 2ND argument must be full."); }
231  it1 = mxGetM(CA_cell_pr);
232  it2 = mxGetN(CA_cell_pr);
233  if (it1 != 1){
234  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
235  mexErrMsgTxt("DSDP: Third column in 2ND argument must have 1 row,");}
236  if (it2 != dim1){
237  mexPrintf("??? DSDP DATA ERROR: Cell: %d \n",j+1);
238  mexErrMsgTxt("DSDP: Secord and third column must have same dimension,");}
239  } else {
240  mexPrintf("??? DSDP DATA ERROR: Cell: %d, Conetype: %s \n",j+1,conetype);
241  mexErrMsgTxt("DSDP: Unknown Cone type in 2ND argument. Try 'SDP' or 'LP' or 'Bounds'. ");
242  }
243  }
244 
245 
246  /* Create output arrays */
247  if (nlhs>2){
248  if (X_OUT != NULL) mxDestroyArray(X_OUT) ;
249  X_OUT = mxCreateCellMatrix(mC,1);
250  }
251  if (Y_OUT != NULL) mxDestroyArray(Y_OUT) ;
252  Y_OUT = mxCreateDoubleMatrix(nvars, 1, mxREAL) ;
253 
254  /* Create the Solver */
255  info = DSDPCreate(nvars,&dsdp); CHKERR(info);
256  info = DSDPCreateSDPCone(dsdp,nsdpblocks,&sdpcone); CHKERR(info);
257 
258  /* Set Dual Objective Vector */
259  bval=mxGetPr(B_IN);
260  if (!bval){ mexErrMsgTxt("DSDP: Problems with 1ST argument");}
261  for (i=0;i<nvars;i++){info=DSDPSetDualObjective(dsdp,i+1,bval[i]); CHKERR(info);}
262 
263  /* Set Matrix Data */
264  for (j=0; j<mC; j++){ /* Begin Block */
265  subs[0] = j; subs[1] = 0;
266  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
267  CA_cell_pr = mxGetCell(CA_IN,index);
268  mxGetString(CA_cell_pr,conetype,20);
269  if (strcmp(conetype,"SDP")==0){
270 
271  subs[0] = j; subs[1] = 1;
272  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
273  CA_cell_pr = mxGetCell(CA_IN,index);
274  it1 = mxGetM(CA_cell_pr);
275  it2 = mxGetN(CA_cell_pr);
276  blockn=mxGetPr(CA_cell_pr);
277  nsubblocks=it1*it2;
278 
279  subs[0] = j; subs[1] = 2;
280  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
281  CA_cell_pr = mxGetCell(CA_IN,index);
282  aval=mxGetPr(CA_cell_pr); air =mxGetIr(CA_cell_pr); ajc =mxGetJc(CA_cell_pr);
283  if (!aval||!air||!ajc)
284  { mexErrMsgTxt("DSDP: Problems with 2ND argument");}
285  for (tnnz=0,jj=0;jj<nsubblocks;jj++){n=(int)blockn[jj];tnnz+=n*(n+1)/2;}
286  if (nlhs>2){
287  subs[0] = j; subs[1] = 0;
288  index = mxCalcSingleSubscript(X_OUT,nsubs,subs);
289  X_cell_pr = mxCreateDoubleMatrix(tnnz,1,mxREAL);
290  mxSetCell(X_OUT,index,X_cell_pr);
291  xout=mxGetPr(X_cell_pr);
292  if (tnnz>0 && !xout){ mexErrMsgTxt("DSDP: Cannot create array. Out of Memory");}
293  }
294 
295  for (ii=0; ii<=nvars; ii++){ /* Begin Variable matrix constraints */
296  i=ii+1;
297  if (i==nvars+1) i=0;
298  tnnz=0; spot=ajc[ii]; blockj=sdpblockj;
299  for (jj=0;jj<nsubblocks;jj++){
300  n=(int)blockn[jj];
301  if (sdpnmax<n) sdpnmax=n;
302  if (ii==0){
303  nn+=n;
304  info=SDPConeSetBlockSize(sdpcone,blockj,n); CHKERR(info);
305  info=SDPConeUsePackedFormat(sdpcone,blockj); CHKERR(info);
306 
307  if (nlhs>2){info=SDPConeSetXArray(sdpcone,blockj,n,xout+tnnz,(n*n+n)/2); CHKERR(info);}
308  info=CountNonzeroMatrices(blockn,nsubblocks,jj,nvars, air, ajc, &nzmats); CHKERR(info);
309  info=SDPConeSetSparsity(sdpcone,blockj,nzmats); CHKERR(info);
310  if (stat1<nzmats)stat1=nzmats;
311  }
312  for (tnnz2=tnnz+n*(n+1)/2,ijnnz=0;ijnnz<ajc[ii+1]-spot && air[spot+ijnnz]<tnnz2;ijnnz++){}
313  if ( ijnnz==0 ){ /* info=DSDPSetZeroMat(dsdp,sdpblockj,i,n); */
314  } else if(ijnnz==n*(n+1)/2){ /* check for dense matrix */
315  if (CheckForConstantMat(aval+spot,ijnnz,n)){
316  info=SDPConeSetConstantMat(sdpcone,blockj,i,n,aval[spot]); CHKERR(info);
317  } else {
318  info=SDPConeSetADenseVecMat(sdpcone,blockj,i,n,1.0,aval+spot,ijnnz); CHKERR(info);
319  }
320  } else {
321  info=SDPConeSetASparseVecMat(sdpcone,blockj,i,n,1.0,tnnz,air+spot,aval+spot,ijnnz); CHKERR(info);
322  }
323  tnnz+=n*(n+1)/2; spot+=ijnnz; blockj++;
324  }
325  } /* End Matrices in block */
326  sdpblockj=sdpblockj+nsubblocks;
327 
328  } else if (strcmp(conetype,"LB")==0 || strcmp(conetype,"UB")==0 ){
329  subs[0] = j; subs[1] = 2;
330  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
331  CA_cell_pr = mxGetCell(CA_IN,index);
332  aval=mxGetPr(CA_cell_pr);
333  info=DSDPCreateBCone(dsdp,&bcone); CHKERR(info);
334  info=BConeAllocateBounds(bcone,nvars); CHKERR(info);
335  for (i=0;i<nvars;i++){
336  if (strcmp(conetype,"LB")==0){
337  info=BConeSetLowerBound(bcone,i+1,aval[i]); CHKERR(info);
338  } else {
339  info=BConeSetUpperBound(bcone,i+1,aval[i]); CHKERR(info);
340  }
341  }
342  if (nlhs>2){
343  subs[0] = j; subs[1] = 0;
344  index = mxCalcSingleSubscript(X_OUT,nsubs,subs);
345  X_cell_pr = mxCreateDoubleMatrix(nvars,1,mxREAL);
346  mxSetCell(X_OUT,index,X_cell_pr);
347  aval2=mxGetPr(X_cell_pr);
348  info=BConeSetXArray(bcone,aval2,nvars); CHKERR(info);
349  }
350  } else if (strcmp(conetype,"LP")==0){
351  subs[0] = j; subs[1] = 2;
352  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
353  CA_cell_pr = mxGetCell(CA_IN,index);
354  n = mxGetM(CA_cell_pr);
355  if (lpnmax<n) lpnmax=n;
356  nn+=n;
357  aval=mxGetPr(CA_cell_pr); air =mxGetIr(CA_cell_pr); ajc =mxGetJc(CA_cell_pr);
358  if (!aval||!air||!ajc)
359  { mexErrMsgTxt("DSDP: Problems with 2ND argument");}
360  info=DSDPCreateLPCone(dsdp,&lpcone); CHKERR(info);
361  info=LPConeSetData2(lpcone,n,ajc,air,aval); CHKERR(info);
362  if (nlhs>2){
363  subs[0] = j; subs[1] = 0;
364  index = mxCalcSingleSubscript(X_OUT,nsubs,subs);
365  X_cell_pr = mxCreateDoubleMatrix(n,1,mxREAL);
366  mxSetCell(X_OUT,index,X_cell_pr);
367  xout=mxGetPr(X_cell_pr);
368  if (n>0 && !xout){ mexErrMsgTxt("DSDP: Cannot create array. Out of Memory");}
369  info=LPConeSetXVec(lpcone,xout,n); CHKERR(info);
370  }
371 
372  } else if (strcmp(conetype,"FIXED")==0){
373  int vari;
374  subs[0] = j; subs[1] = 1;
375  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
376  CA_cell_pr = mxGetCell(CA_IN,index);
377  aval=mxGetPr(CA_cell_pr);
378  it1 = mxGetM(CA_cell_pr);
379  it2 = mxGetN(CA_cell_pr);
380  subs[0] = j; subs[1] = 2;
381  index = mxCalcSingleSubscript(CA_IN,nsubs,subs);
382  CA_cell_pr = mxGetCell(CA_IN,index);
383  aval2=mxGetPr(CA_cell_pr);
384  if (nlhs>2){
385  subs[0] = j; subs[1] = 0;
386  index = mxCalcSingleSubscript(X_OUT,nsubs,subs);
387  X_cell_pr = mxCreateDoubleMatrix(it1*it2,1,mxREAL);
388  mxSetCell(X_OUT,index,X_cell_pr);
389  xout=mxGetPr(X_cell_pr);
390  if (it1*it2>0 && !xout){ mexErrMsgTxt("DSDP: Cannot create array. Out of Memory");}
391  } else {xout=0;}
392  info=DSDPSetFixedVariables(dsdp,aval,aval2,xout,it1*it2); CHKERR(info);
393  for (i=0;i<it1*it2;i++){
394  /*
395  vari=(int)aval[i];
396  printf("FixedVariable %d to %4.4e\n",vari,aval2[i]);
397  info=DSDPSetFixedVariable(dsdp,vari,aval2[i]); CHKERR(info);
398  */
399  }
400  }
401  } /* End Block */
402 
403  /* Set initial point */
404  if (y0){
405  for (i=0;i<nvars;i++){ info = DSDPSetY0(dsdp,i+1,y0[i]); CHKERR(info);}
406  }
407 
408  reuse=(nvars-2)/sdpnmax;
409  if (nvars<50 && reuse==0) reuse=1;
410  if (reuse>=1) reuse++;
411  reuse=reuse*reuse;
412  if (nvars<2000 && reuse>10) reuse=10;
413  if (reuse>12) reuse=12;
414  info=DSDPReuseMatrix(dsdp,reuse); CHKERR(info);
415 
416  info=DSDPGetPPObjective(dsdp,&zbar); CHKERR(info);
417  info=DSDPGetR(dsdp,&r0); CHKERR(info);
418  info=DSDPGetMaxIts(dsdp,&maxit); CHKERR(info);
419  info=DSDPGetPenaltyParameter(dsdp,&penalty); CHKERR(info);
420  info=DSDPGetPotentialParameter(dsdp,&rho); CHKERR(info);
421  info=DSDPGetDualBound(dsdp,&dbound); CHKERR(info);
422  info=DSDPGetGapTolerance(dsdp,&gaptol); CHKERR(info);
423  info=DSDPGetRTolerance(dsdp,&inftol); CHKERR(info);
424  info=DSDPGetBarrierParameter(dsdp,&mu0); CHKERR(info);
425  info=DSDPGetMaxTrustRadius(dsdp,&maxtrust); CHKERR(info);
426  info=DSDPGetStepTolerance(dsdp,&steptol); CHKERR(info);
427  info=DSDPGetPTolerance(dsdp,&infptol); CHKERR(info);
428  info=DSDPGetDataNorms(dsdp, datanorm); CHKERR(info);
429  info=DSDPGetYBounds(dsdp,&ylow,&yhigh); CHKERR(info);
430  if (datanorm[0]==0){info=DSDPSetYBounds(dsdp,-1.0,1.0);CHKERR(info);}
431  if (nrhs >2){
432  if(!mxIsStruct(PARS_IN)){
433  mexErrMsgTxt("Fifth Parameter `OPTIONS' should be a structure.");}
434 
435  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"maxit") ){
436  maxit=(int) mxGetScalar(OPTIONS_FIELD);
437  info=DSDPSetMaxIts(dsdp,maxit); CHKERR(info);}
438  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"fastblas") ){
439  fastblas= (int) mxGetScalar(OPTIONS_FIELD);}
440  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"print") ){
441  print_info= (int) mxGetScalar(OPTIONS_FIELD);
442  printlevel=print_info;
443  }
444  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"bigM") ){
445  rpos=(int)mxGetScalar(OPTIONS_FIELD);
446  info=DSDPUsePenalty(dsdp,rpos); CHKERR(info);}
447  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"penalty") ){
448  penalty=mxGetScalar(OPTIONS_FIELD);
449  info=DSDPSetPenaltyParameter(dsdp,penalty); CHKERR(info);}
450  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"rho") ){
451  rho= mxGetScalar(OPTIONS_FIELD);
452  info=DSDPSetPotentialParameter(dsdp,rho); CHKERR(info);}
453  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"dynamicrho") ){
454  drho= (int)mxGetScalar(OPTIONS_FIELD);
455  info=DSDPUseDynamicRho(dsdp,drho); CHKERR(info);}
456  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"zbar") ){
457  zbar= mxGetScalar(OPTIONS_FIELD);
458  info=DSDPSetZBar(dsdp,zbar); CHKERR(info);}
459  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"dual_bound") ){
460  dbound= mxGetScalar(OPTIONS_FIELD);
461  info=DSDPSetDualBound(dsdp,dbound); CHKERR(info);}
462  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"reuse") ){
463  reuse= (int)mxGetScalar(OPTIONS_FIELD);
464  info=DSDPReuseMatrix(dsdp,reuse); CHKERR(info);}
465  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"gaptol") ){
466  gaptol= mxGetScalar(OPTIONS_FIELD);
467  info=DSDPSetGapTolerance(dsdp,gaptol);CHKERR(info);}
468  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"lp_barrier") ){
469  lpb= mxGetScalar(OPTIONS_FIELD);
470  if (lpb<0.1) lpb=0.1;
471  if (lpcone){info=LPConeScaleBarrier(lpcone,lpb); CHKERR(info);}}
472  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"lpb") ){
473  lpb= mxGetScalar(OPTIONS_FIELD);
474  if (lpb<0.1) lpb=0.1;
475  if (lpcone){info=LPConeScaleBarrier(lpcone,lpb); CHKERR(info);}}
476  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"cc") ){
477  cc= mxGetScalar(OPTIONS_FIELD);
478  info=DSDPAddObjectiveConstant(dsdp,cc); CHKERR(info);}
479  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"inftol") ){
480  inftol= mxGetScalar(OPTIONS_FIELD);
481  info=DSDPSetRTolerance(dsdp,inftol); CHKERR(info);}
482  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"infptol") ){
483  infptol= mxGetScalar(OPTIONS_FIELD);
484  info=DSDPSetPTolerance(dsdp,infptol);CHKERR(info);}
485  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"pnormtol") ){
486  pnormtol= mxGetScalar(OPTIONS_FIELD);
487  info=DSDPSetPNormTolerance(dsdp,pnormtol);CHKERR(info);}
488  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"boundy") ){
489  yhigh= fabs(mxGetScalar(OPTIONS_FIELD)); ylow=-yhigh;
490  info=DSDPSetYBounds(dsdp,ylow,yhigh);CHKERR(info);}
491  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"r0") ){
492  r0= mxGetScalar(OPTIONS_FIELD);
493  info=DSDPSetR0(dsdp,r0); CHKERR(info);}
494  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"mu0") ){
495  mu0=mxGetScalar(OPTIONS_FIELD);
496  info=DSDPSetBarrierParameter(dsdp,mu0);CHKERR(info);}
497  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"maxtrustradius")){
498  maxtrust= mxGetScalar(OPTIONS_FIELD);
499  info=DSDPSetMaxTrustRadius(dsdp,maxtrust); CHKERR(info);}
500  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"steptol") ){
501  steptol = mxGetScalar(OPTIONS_FIELD);
502  info=DSDPSetStepTolerance(dsdp,steptol); CHKERR(info);}
503  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"dobjmin") ){
504  dtmp= mxGetScalar(OPTIONS_FIELD);
505  info = DSDPSetDualLowerBound(dsdp,dtmp);}
506  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"dloginfo") ){
507  iloginfo= (int) mxGetScalar(OPTIONS_FIELD);
508  info=DSDPLogInfoAllow(iloginfo,0);}
509  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"logtime") ){
510  printsummary= (int) mxGetScalar(OPTIONS_FIELD);}
511  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"printproblem") ){
512  info=DSDPPrintData(dsdp,sdpcone,lpcone);}
513  if ( OPTIONS_FIELD = mxGetField(PARS_IN,0,"xmaker") ){
514  xmaker = (int) mxGetScalar(OPTIONS_FIELD);}
515  /*
516  if( OPTIONS_FIELD = mxGetField(PARS_IN,0,"maxlanczos") ){
517  itmp= (int) mxGetScalar(OPTIONS_FIELD);
518  info=DSDPSetLanczosIterations(dsdp,itmp);}
519  */
520 
521  }
522 
523  info = DSDPSetMonitor(dsdp,DSDPPrintStats2,0); CHKERR(info);
524 
525  info = DSDPSetup(dsdp);
526  if (info){
527  mexErrMsgTxt("DSDP: Setup Error, Probably out of memory");}
528 
529 
530  info = DSDPSolve(dsdp); CHKERR(info);
531  info = DSDPStopReason(dsdp,&reason);
532  if (reason!=DSDP_INFEASIBLE_START){
533  info=DSDPComputeX(dsdp);CHKERR(info);
534  }
535  info = DSDPGetSolutionType(dsdp,&pdfeasible); CHKERR(info);
536 
537  if (printsummary){ DSDPEventLogSummary();}
538 
539  if (info){
540  mexErrMsgTxt("DSDP: Numerical error");}
541 
542  if ( reason == DSDP_INFEASIBLE_START){
543  mexErrMsgTxt("DSDP Terminated Due to Infeasible Starting Point\n");
544  } else if (print_info){
545 
546  if (reason == DSDP_CONVERGED)
547  mexPrintf("DSDP Converged. \n");
548  else if ( reason == DSDP_UPPERBOUND )
549  mexPrintf("DSDP Converged: Dual Objective exceeds its bound\n");
550  else if ( reason == DSDP_SMALL_STEPS )
551  mexPrintf("DSDP Terminated Due to Small Steps\n");
552  else if ( reason == DSDP_MAX_IT)
553  mexPrintf("DSDP Terminated Due Maximum Number of Iterations\n");
554  else if ( reason == DSDP_USER_TERMINATION)
555  mexPrintf("DSDP Terminated By User\n");
556  else if ( reason == DSDP_INFEASIBLE_START)
557  mexPrintf("DSDP Terminated Due to Infeasible Starting Point\n");
558  else
559  mexPrintf("DSDP Finished.\n");
560  }
561 
562  if (pdfeasible==DSDP_UNBOUNDED){
563  mexPrintf("DSDP: Dual Unbounded, Primal Infeasible\n");
564  } else if (pdfeasible==DSDP_INFEASIBLE){
565  mexPrintf("DSDP: Primal Unbounded, Dual Infeasible\n");
566  }
567 
568  /* Set the dual solution */
569  yout=mxGetPr(Y_OUT);
570  info = DSDPGetY(dsdp,yout,nvars); CHKERR(info);
571  if (info){
572  mexErrMsgTxt("DSDP: Numerical error");}
573 
574 
575  /* Output statistics */
576  if (STAT_OUT != NULL) mxDestroyArray(STAT_OUT) ;
577  subs[0] = 1; subs[1] = 1;
578  STAT_OUT = mxCreateStructArray(2,subs,nfields,fnames);
579  info= DSDPGetDObjective(dsdp,&dobj); CHKERR(info);
580  info= DSDPGetPObjective(dsdp,&pobj); CHKERR(info);
581  info= DSDPGetR(dsdp,&dinf); CHKERR(info);
582  info= DSDPStopReason(dsdp,&reason); CHKERR(info);
583  info= DSDPGetIts(dsdp,&its); CHKERR(info);
584  info= DSDPGetBarrierParameter(dsdp,&mu); CHKERR(info);
585  info= DSDPGetStepLengths(dsdp,&pstep,&dstep); CHKERR(info);
586  info= DSDPGetPnorm(dsdp,&pnorm); CHKERR(info);
587  info= DSDPGetBarrierParameter(dsdp,&mu); CHKERR(info);
588  info= DSDPGetYBounds(dsdp,&ylow,&yhigh); CHKERR(info);
589 
590  if (pdfeasible==DSDP_UNBOUNDED){
591  STAT_FIELD = mxCreateString("Unbounded");
592  mxSetField(STAT_OUT,0,fnames[0],STAT_FIELD);
593  } else if (pdfeasible==DSDP_INFEASIBLE){
594  STAT_FIELD = mxCreateString("Infeasible");
595  mxSetField(STAT_OUT,0,fnames[0],STAT_FIELD);
596  } else {
597  STAT_FIELD = mxCreateString("PDFeasible");
598  mxSetField(STAT_OUT,0,fnames[0],STAT_FIELD);
599  }
600 
601  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
602  stat=mxGetPr(STAT_FIELD); stat[0]=dobj;
603  mxSetField(STAT_OUT,0,fnames[1],STAT_FIELD);
604 
605  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
606  stat=mxGetPr(STAT_FIELD); stat[0]=pobj;
607  mxSetField(STAT_OUT,0,fnames[2],STAT_FIELD);
608 
609  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
610  stat=mxGetPr(STAT_FIELD); stat[0]=dobj;
611  mxSetField(STAT_OUT,0,fnames[3],STAT_FIELD);
612 
613  if (reason==DSDP_CONVERGED){
614  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
615  stat=mxGetPr(STAT_FIELD); stat[0]=0;
616  mxSetField(STAT_OUT,0,fnames[4],STAT_FIELD);
617  } else {
618  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
619  stat=mxGetPr(STAT_FIELD); stat[0]=(double)reason;
620  if (stat[0]==0) stat[0]=-1;
621  mxSetField(STAT_OUT,0,fnames[4],STAT_FIELD);
622  }
623 
624  if (reason==DSDP_CONVERGED){ /* For YALMIP */
625  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
626  stat=mxGetPr(STAT_FIELD); stat[0]=0;
627  if (pdfeasible==DSDP_UNBOUNDED){ stat[0]=1;}
628  if (pdfeasible==DSDP_INFEASIBLE){ stat[0]=2;}
629  mxSetField(STAT_OUT,0,fnames[5],STAT_FIELD);
630  } else {
631  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
632  stat=mxGetPr(STAT_FIELD); stat[0]=-1;
633  if (reason==DSDP_INFEASIBLE_START)stat[0]=-6;
634  if (reason==DSDP_UPPERBOUND)stat[0]=-5;
635  if (reason==DSDP_MAX_IT)stat[0]=-3;
636  mxSetField(STAT_OUT,0,fnames[5],STAT_FIELD);
637  }
638 
639  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
640  stat=mxGetPr(STAT_FIELD); stat[0]=(double) its;
641  mxSetField(STAT_OUT,0,fnames[6],STAT_FIELD);
642 
643  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
644  stat=mxGetPr(STAT_FIELD); stat[0]=dinf;
645  mxSetField(STAT_OUT,0,fnames[7],STAT_FIELD);
646 
647  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
648  stat=mxGetPr(STAT_FIELD); stat[0]=mu;
649  mxSetField(STAT_OUT,0,fnames[8],STAT_FIELD);
650 
651  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
652  stat=mxGetPr(STAT_FIELD); stat[0]=pstep;
653  mxSetField(STAT_OUT,0,fnames[9],STAT_FIELD);
654 
655  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
656  stat=mxGetPr(STAT_FIELD); stat[0]=dstep;
657  mxSetField(STAT_OUT,0,fnames[10],STAT_FIELD);
658 
659  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
660  stat=mxGetPr(STAT_FIELD); stat[0]=pnorm;
661  mxSetField(STAT_OUT,0,fnames[11],STAT_FIELD);
662 
663  itmp=100; if (its < itmp) itmp=its;
664  STAT_FIELD = mxCreateDoubleMatrix(itmp, 1, mxREAL) ;
665  stat=mxGetPr(STAT_FIELD);
666  info= DSDPGetGapHistory(dsdp,stat,itmp); CHKERR(info);
667  mxSetField(STAT_OUT,0,fnames[12],STAT_FIELD);
668 
669  STAT_FIELD = mxCreateDoubleMatrix(itmp, 1, mxREAL) ;
670  stat=mxGetPr(STAT_FIELD);
671  info= DSDPGetRHistory(dsdp,stat,itmp); CHKERR(info);
672  mxSetField(STAT_OUT,0,fnames[13],STAT_FIELD);
673 
674  STAT_FIELD = mxCreateDoubleMatrix(6, 1, mxREAL) ;
675  stat=mxGetPr(STAT_FIELD);
676  info = DSDPGetFinalErrors(dsdp,stat); CHKERR(info);
677  mxSetField(STAT_OUT,0,fnames[14],STAT_FIELD);
678 
679  STAT_FIELD = mxCreateDoubleMatrix(3, 1, mxREAL) ;
680  stat=mxGetPr(STAT_FIELD);
681  info = DSDPGetDataNorms(dsdp,stat); CHKERR(info);
682  dtmp=stat[0];stat[0]=stat[1];stat[1]=dtmp;
683  mxSetField(STAT_OUT,0,fnames[15],STAT_FIELD);
684 
685  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
686  stat=mxGetPr(STAT_FIELD);
687  info = DSDPGetYMaxNorm(dsdp,stat); CHKERR(info);
688  mxSetField(STAT_OUT,0,fnames[16],STAT_FIELD);
689 
690  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
691  stat=mxGetPr(STAT_FIELD);
692  stat[0]=yhigh;
693  mxSetField(STAT_OUT,0,fnames[17],STAT_FIELD);
694 
695  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
696  stat=mxGetPr(STAT_FIELD);
697  info = DSDPGetPenaltyParameter(dsdp,stat); CHKERR(info);
698  mxSetField(STAT_OUT,0,fnames[18],STAT_FIELD);
699 
700  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
701  stat=mxGetPr(STAT_FIELD);
702  info = DSDPGetTraceX(dsdp,stat); CHKERR(info);
703  mxSetField(STAT_OUT,0,fnames[19],STAT_FIELD);
704 
705  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
706  info = DSDPGetReuseMatrix(dsdp,&its); CHKERR(info);
707  stat=mxGetPr(STAT_FIELD); stat[0]=(double) its;
708  mxSetField(STAT_OUT,0,fnames[20],STAT_FIELD);
709 
710  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
711  stat=mxGetPr(STAT_FIELD);
712  info = DSDPGetPotentialParameter(dsdp,stat); CHKERR(info);
713  mxSetField(STAT_OUT,0,fnames[21],STAT_FIELD);
714 
715  if (xmaker){
716  STAT_FIELD = mxCreateDoubleMatrix(nvars+1, 1, mxREAL) ;
717  stat=mxGetPr(STAT_FIELD);
718  info = DSDPGetYMakeX(dsdp,stat,nvars+1); CHKERR(info);
719  mxSetField(STAT_OUT,0,fnames[22],STAT_FIELD);
720 
721  STAT_FIELD = mxCreateDoubleMatrix(nvars+1, 1, mxREAL) ;
722  stat=mxGetPr(STAT_FIELD);
723  info = DSDPGetDYMakeX(dsdp,stat,nvars+1); CHKERR(info);
724  mxSetField(STAT_OUT,0,fnames[23],STAT_FIELD);
725 
726  STAT_FIELD = mxCreateDoubleMatrix(1, 1, mxREAL) ;
727  stat=mxGetPr(STAT_FIELD);
728  info = DSDPGetMuMakeX(dsdp,stat); CHKERR(info);
729  mxSetField(STAT_OUT,0,fnames[24],STAT_FIELD);
730  }
731  /* Free internal data structure */
732 
733  info = DSDPDestroy(dsdp); CHKERR(info);
734 
735  return;
736 } /* main */
737 
738 
739 
740 #undef __FUNCT__
741 #define __FUNCT__ "CheckForConstantMat"
742 static int CheckForConstantMat(double*v,int nnz, int n){
743  int i;double vv;
744  if (n<=1){ return 0; }
745  if (nnz!=(n*n+n)/2){ return 0; }
746  for (vv=v[0],i=1;i<nnz;i++){
747  if (v[i]!=vv){ return 0;}
748  }
749  return 1;
750 }
751 
752 #undef __FUNCT__
753 #define __FUNCT__ "CountNonzeroMatrices"
754  static int CountNonzeroMatrices(double *blocksize, int nblocks, int block, int nvars, int *indd, int*nnnz, int *nnzmats){
755  int i,j,n,nzmats=0;
756  int marker1=0,marker2=0;
757  for (i=0;i<block;i++){n=(int)blocksize[i]; marker1=marker1+n*(n+1)/2;}
758  n=(int)blocksize[block];
759  marker2=marker1+n*(n+1)/2;
760  for (i=0;i<nvars;i++){
761  j=nnnz[i];
762  while (indd[j]<marker1 && j< nnnz[i+1]){j++;}
763  if (j<nnnz[i+1] && indd[j]<marker2){ nzmats++;}
764  }
765  *nnzmats=nzmats;
766  return 0;
767 }
768 
769 /* ---------------------------------------------------------- */
770  static int DSDPPrintStats2(DSDP dsdp, void* dummy){
771 
772  int iter,info;
773  double pobj,dobj,pstp=0,dstp,mu,res,pnorm,pinfeas;
774  DSDPTerminationReason reason;
775 
776  if(printlevel<=0) return(0);
777 
778  info = DSDPStopReason(dsdp,&reason);
779  info = DSDPGetIts(dsdp,&iter);
780 
781  if( (reason!=CONTINUE_ITERATING) || ((iter % printlevel)==0)){
782 
783  info = DSDPGetDDObjective(dsdp,&dobj);
784  info = DSDPGetPPObjective(dsdp,&pobj);
785  info = DSDPGetR(dsdp,&res);
786  info = DSDPGetPInfeasibility(dsdp,&pinfeas);
787  info = DSDPGetStepLengths(dsdp,&pstp,&dstp);
788  info = DSDPGetBarrierParameter(dsdp,&mu);
789  info = DSDPGetPnorm(dsdp,&pnorm);
790 
791 
792  if (iter==0){
793  mexPrintf("Iter PP Objective DD Objective PInfeas DInfeas Nu StepLength Pnrm\n")
794  ;
795  mexPrintf("----------------------------------------------------------------------------------------\n")
796  ;
797  }
798  mexPrintf("%-3d %16.8e %16.8e %9.1e %9.1e %9.1e",iter,pobj,dobj,pinfeas,res,mu);
799  mexPrintf(" %4.2f %4.2f",pstp,dstp);
800  if (pnorm>1.0e3){
801  mexPrintf(" %1.0e \n",pnorm);
802  } else {
803  mexPrintf(" %5.2f \n",pnorm);
804  }
805  fflush(NULL);
806  }
807  return(0);
808 }
int DSDPSetBarrierParameter(DSDP dsdp, double mu)
Set the current barrier parameter.
Definition: dsdpsetdata.c:340
int DSDPCreate(int m, DSDP *dsdpnew)
Create a DSDP solver. FIRST DSDP routine!
Definition: dsdpsetup.c:30
int DSDPGetFinalErrors(DSDP dsdp, double err[6])
Copy six different error measurements into an array.
Definition: dsdpx.c:297
int DSDPGetDualBound(DSDP dsdp, double *dbound)
Get the termination parameter.
Definition: dsdpconverge.c:227
int DSDPDestroy(DSDP dsdp)
Free the internal data structures of the solver and the cones associated with it.
Definition: dsdpsetup.c:496
int DSDPGetStepTolerance(DSDP dsdp, double *steptol)
Get the current tolerance.
Definition: dsdpconverge.c:273
int DSDPGetStepLengths(DSDP dsdp, double *pstep, double *dstep)
Copy the step sizes in the current iteration.
Definition: dsdpsetdata.c:742
int DSDPGetBarrierParameter(DSDP dsdp, double *mu)
Copy the current barrier parameter.
Definition: dsdpsetdata.c:364
int DSDPGetDYMakeX(DSDP dsdp, double dy[], int m)
Copies the variables dy used to construct X into an array.
Definition: dsdpx.c:483
int DSDPGetPPObjective(DSDP dsdp, double *ppobj)
Copy the objective value (PP).
Definition: dsdpsetdata.c:479
int DSDPAddObjectiveConstant(DSDP dsdp, double c)
Add a constant to the objective.
Definition: dsdpsetdata.c:185
int DSDPGetRHistory(DSDP dsdp, double hist[], int length)
Copy a history of the infeasibility in (D) into an array.
Definition: dsdpconverge.c:298
int DSDPSetDualBound(DSDP dsdp, double dbound)
Terminate the solver if the objective value in (DD) is greater than this tolerance.
Definition: dsdpconverge.c:205
int DSDPGetPenaltyParameter(DSDP dsdp, double *Gamma)
Copy the penalty parameter Gamma.
Definition: dsdpsetdata.c:443
int DSDPCreateBCone(DSDP dsdp, BCone *dspcone)
Create a new cone that represents bounds on the y variables.
Definition: dbounds.c:467
int DSDPGetGapTolerance(DSDP dsdp, double *gaptol)
Get the termination tolerance.
Definition: dsdpconverge.c:132
int BConeSetUpperBound(BCone bcone, int vari, double ubound)
Set an upper bound on a variable y.
Definition: dbounds.c:583
int DSDPSetZBar(DSDP dsdp, double ppobj)
Set an upper bound on the objective value at the solution.
Definition: dsdpsetdata.c:283
int DSDPGetPTolerance(DSDP dsdp, double *inftol)
Copy the feasibility tolerance.
Definition: dsdpx.c:386
int DSDPSetDualObjective(DSDP dsdp, int i, double bi)
Set the objective vector b in (D).
Definition: dsdpsetdata.c:25
int DSDPUseDynamicRho(DSDP dsdp, int yesorno)
Use a dynamic strategy to choose parameter rho.
Definition: dsdpsetdata.c:821
int DSDPGetPnorm(DSDP dsdp, double *pnorm)
Copy the proximity of the solution to the central path.
Definition: dsdpsetdata.c:724
int SDPConeSetXArray(SDPCone sdpcone, int blockj, int n, double xx[], int nn)
Provide an array for the SDPCone object can use to store dense matrices.
Definition: dsdpadddata.c:278
Internal structures for the DSDP solver.
Definition: dsdp.h:65
int DSDPCreateLPCone(DSDP dsdp, LPCone *dspcone)
Create a new object for linear programs and scalar inequalities.
Definition: dsdplp.c:509
DSDPTerminationReason
There are many reasons to terminate the solver.
int DSDPSetGapTolerance(DSDP dsdp, double gaptol)
Terminate the solver when the relative duality gap is less than this tolerance.
Definition: dsdpconverge.c:110
int DSDPPrintData(DSDP dsdp, SDPCone sdpcone, LPCone lpcone)
Print data in SDPA format to a file named "output.sdpa".
Definition: printsdpa.c:164
int DSDPGetDataNorms(DSDP dsdp, double dnorm[3])
Copy the norms of the data C, A, and b into an array.
Definition: dsdpsetdata.c:621
int BConeSetLowerBound(BCone bcone, int vari, double lbound)
Set a lower bound on a variable y.
Definition: dbounds.c:566
The API to DSDP for those applications using DSDP as a subroutine library.
int DSDPSetRTolerance(DSDP dsdp, double inftol)
Classify (D) as feasible only if the variable r is less than this tolerance.
Definition: dsdpx.c:409
int DSDPSetPTolerance(DSDP dsdp, double inftol)
Classify (P) as feasible only if the infeasibility is less than this tolerance.
Definition: dsdpx.c:365
int SDPConeSetConstantMat(SDPCone sdpcone, int blockj, int vari, int n, double value)
Set a matrix whose elements are all the same.
int DSDPSetMaxTrustRadius(DSDP dsdp, double rad)
Set a maximum trust radius on the step direction.
Definition: dsdpsetdata.c:246
int DSDPGetSolutionType(DSDP dsdp, DSDPSolutionType *pdfeasible)
Solutions can be bounded, infeasible, or unbounded.
Definition: dsdpx.c:254
int DSDPSetup(DSDP dsdp)
Set up data structures in the solver and the cones associated with it.
Definition: dsdpsetup.c:193
int DSDPReuseMatrix(DSDP dsdp, int rm)
Reuse the Hessian of the barrier function multiple times at each DSDP iteration.
Definition: dsdpsetdata.c:905
int BConeAllocateBounds(BCone bcone, int nnz)
Set a surplus variable in constraint in (P).
Definition: dbounds.c:645
int DSDPSetPenaltyParameter(DSDP dsdp, double Gamma)
Set the penalty parameter Gamma.
Definition: dsdpsetdata.c:418
int DSDPGetReuseMatrix(DSDP dsdp, int *rm)
Copy this parameter.
Definition: dsdpsetdata.c:925
int DSDPGetPInfeasibility(DSDP dsdp, double *pperror)
Copy the infeasibility in (P).
Definition: dsdpx.c:343
int DSDPGetMaxTrustRadius(DSDP dsdp, double *rad)
Copy the current radius of the trust region.
Definition: dsdpsetdata.c:265
int DSDPGetDObjective(DSDP dsdp, double *dobj)
Copy the objective value (D).
Definition: dsdpsetdata.c:502
int DSDPSetFixedVariables(DSDP dsdp, double vars[], double vals[], double xout[], int nvars)
Fix variable y to exact values.
Definition: dsdpschurmat.c:695
int DSDPSetPotentialParameter(DSDP dsdp, double rho)
Set the potential parameter.
Definition: dsdpsetdata.c:765
int DSDPStopReason(DSDP dsdp, DSDPTerminationReason *reason)
Copy the reason why the solver terminated.
Definition: dsdpsetdata.c:582
int DSDPGetPObjective(DSDP dsdp, double *pobj)
Copy the objective value (P).
Definition: dsdpx.c:232
int DSDPSetY0(DSDP dsdp, int i, double yi0)
Set the initial values of variables y in (D).
Definition: dsdpsetdata.c:77
int DSDPSetStepTolerance(DSDP dsdp, double steptol)
Terminate the solver if the step length in (DD) is below this tolerance.
Definition: dsdpconverge.c:252
int SDPConeUsePackedFormat(SDPCone sdpcone, int blockj)
Use packed symmetric format for the dense array.
Definition: dsdpadddata.c:452
int DSDPGetIts(DSDP dsdp, int *its)
Copy the current iteration number.
Definition: dsdpsetdata.c:564
int DSDPSetPNormTolerance(DSDP dsdp, double ptol)
Terminate the solver when the relative duality gap is suffiently small and the PNorm is less than thi...
Definition: dsdpconverge.c:158
int DSDPGetMuMakeX(DSDP dsdp, double *mu)
Copies the value of mu used to construct X.
Definition: dsdpx.c:511
int DSDPGetYMaxNorm(DSDP dsdp, double *ynorm)
Copy the the infinity norm of the variables y.
Definition: dsdpsetdata.c:678
int DSDPSolve(DSDP dsdp)
Apply DSDP to the problem.
Definition: dsdpsetup.c:343
int SDPConeSetBlockSize(SDPCone sdpcone, int blockj, int n)
Set the dimension of one block in the semidefinite cone.
Definition: dsdpadddata.c:535
Internal structure for semidefinite cone.
Definition: dsdpsdp.h:80
int DSDPGetY(DSDP dsdp, double y[], int m)
Copies the variables y into an array.
Definition: dsdpsetdata.c:100
int SDPConeSetSparsity(SDPCone sdpcone, int blockj, int nnz)
Set the number of nonzero matrices in a block of the semidefinite cone.
Definition: dsdpadddata.c:596
int DSDPUsePenalty(DSDP dsdp, int yesorno)
Use penalty parameter to enforce feasibility.
Definition: dsdpsetdata.c:383
int DSDPGetDDObjective(DSDP dsdp, double *ddobj)
Copy the objective value (DD).
Definition: dsdpsetdata.c:523
int DSDPGetRTolerance(DSDP dsdp, double *inftol)
Copy the maximum infeasibility allowed (D).
Definition: dsdpx.c:434
int DSDPGetTraceX(DSDP dsdp, double *tracex)
Copy the trace of the variables X in (P).
Definition: dsdpx.c:278
int DSDPGetR(DSDP dsdp, double *res)
Copy the infeasibility in (D), or the variable r in (DD).
Definition: dsdpsetdata.c:601
int DSDPComputeX(DSDP dsdp)
Compute the X variables.
Definition: dsdpx.c:55
int DSDPGetMaxIts(DSDP dsdp, int *its)
Copy the maximum number of iterations from the solver.
Definition: dsdpsetdata.c:225
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
Call DSDP from the Matlab environment.
Definition: dsdp.c:35
int SDPConeSetASparseVecMat(SDPCone sdpcone, int blockj, int vari, int n, double alpha, int ishift, const int ind[], const double val[], int nnz)
Set data matrix in a sparse format.
struct BCone_C * BCone
The BCone object points to lower and upper bounds on the variable y in (D).
Definition: dsdp5.h:28
int SDPConeSetADenseVecMat(SDPCone sdpcone, int blockj, int vari, int n, double alpha, double val[], int nnz)
Set a matrix in a dense format.
int DSDPSetYBounds(DSDP dsdp, double lbound, double ubound)
Bound the variables y.
Definition: dsdpsetdata.c:885
int DSDPGetPotentialParameter(DSDP dsdp, double *rho)
Copy the potential parameter.
Definition: dsdpsetdata.c:784
int DSDPGetYBounds(DSDP dsdp, double *lbound, double *ubound)
Copy the bounds on the variables y.
Definition: dsdpsetdata.c:866
int LPConeSetData2(LPCone lpcone, int n, const int ik[], const int cols[], const double vals[])
Set data A and into the LP cone.
Definition: dsdplp.c:717
int DSDPSetR0(DSDP dsdp, double res)
Set an initial value for the variable r in (DD)
Definition: dsdpsetdata.c:311
int DSDPSetMaxIts(DSDP dsdp, int its)
Terminate the solver after this number of iterations.
Definition: dsdpsetdata.c:206
int DSDPGetGapHistory(DSDP dsdp, double hist[], int length)
Copy a history of the duality gap into an array.
Definition: dsdpconverge.c:321
int DSDPGetYMakeX(DSDP dsdp, double y[], int m)
Copies the variables y used to construct X into an array.
Definition: dsdpx.c:455
struct LPCone_C * LPCone
The LPCone object points to blocks of data that specify linear scalar inequality constraints.
Definition: dsdp5.h:27
DSDPSolutionType
Formulations (P) and (D) can be feasible and bounded, feasible and unbounded, or infeasible.