DSDP
dsdpsetup.c
Go to the documentation of this file.
1#include "dsdp.h"
2#include "dsdpsys.h"
3#include "dsdp5.h"
8
28#undef __FUNCT__
29#define __FUNCT__ "DSDPCreate"
30int DSDPCreate(int m,DSDP* dsdpnew){
31
32 DSDP dsdp;
33 int info;
34
35 DSDPFunctionBegin;
36
37 DSDPCALLOC1(&dsdp,PD_DSDP,&info);DSDPCHKERR(info);
38 *dsdpnew=dsdp;
39 dsdp->keyid=DSDPKEY;
40
41 /* Initialize some parameters */
42 DSDPEventLogInitialize();
43 dsdp->m=m;
44 dsdp->maxcones=0;
45 dsdp->ncones=0;
46 dsdp->K=0;
47 dsdp->setupcalled=DSDP_FALSE;
48 dsdp->ybcone=0;
49 dsdp->ndroutines=0;
50 /* info = DSDPSetStandardMonitor(dsdp);DSDPCHKERR(info); */
51 info = DSDPVecCreateSeq(m+2,&dsdp->b);DSDPCHKERR(info);
52 info = DSDPVecZero(dsdp->b);DSDPCHKERR(info);
53 info = DSDPVecDuplicate(dsdp->b,&dsdp->y);DSDPCHKERR(info);
54 info = DSDPVecDuplicate(dsdp->b,&dsdp->ytemp);DSDPCHKERR(info);
55 info = DSDPVecZero(dsdp->y);DSDPCHKERR(info);
56 info = DSDPVecSetC(dsdp->y,-1.0);DSDPCHKERR(info);
57
58 info = DSDPAddRCone(dsdp,&dsdp->rcone);DSDPCHKERR(info);
59 info = DSDPCreateLUBoundsCone(dsdp,&dsdp->ybcone);DSDPCHKERR(info);
60
61 info=DSDPSetDefaultStatistics(dsdp);DSDPCHKERR(info);
62 info=DSDPSetDefaultParameters(dsdp);DSDPCHKERR(info);
63 info=DSDPSetDefaultMonitors(dsdp);DSDPCHKERR(info);
64
65 /* info = DSDPMatInitialize(m,m,&dsdp->Q);DSDPCHKERR(info); */
66 info = DSDPSchurMatInitialize(&dsdp->M);DSDPCHKERR(info);
67 info = DSDPSetDefaultSchurMatrixStructure(dsdp); DSDPCHKERR(info);
68 info = DSDPCGInitialize(&dsdp->sles); DSDPCHKERR(info);
69
70 /* Set the one global variable
71 sdat=dsdp;
72 */
73 DSDPFunctionReturn(0);
74}
75
76
77#undef __FUNCT__
78#define __FUNCT__ "DSDPSetDefaultStatistics"
85
86 int i;
87 DSDPFunctionBegin;
88 DSDPValid(dsdp);
89 dsdp->reason=CONTINUE_ITERATING;
90 dsdp->pdfeasible=DSDP_PDUNKNOWN;
91 dsdp->itnow=0;
92 dsdp->pobj= 1.0e10;
93 dsdp->ppobj= 1.0e10;
94 dsdp->dobj= -1.0e+9;
95 dsdp->ddobj= -1.0e+9;
96 dsdp->dualitygap=dsdp->ppobj-dsdp->ddobj;
97 dsdp->pstep=1.0;
98 dsdp->dstep=0.0;
99 for (i=0;i<MAX_XMAKERS;i++){
100 dsdp->xmaker[i].mu=1.0e200;
101 dsdp->xmaker[i].pstep=0.0;
102 }
103 dsdp->pnorm=0.001;
104 dsdp->mu=1000.0;
105 dsdp->np=0;
106 dsdp->anorm=0;
107 dsdp->bnorm=0;
108 dsdp->cnorm=0;
109 dsdp->tracex=0;
110 dsdp->tracexs=0;
111 dsdp->Mshift=0;
112 dsdp->goty0=DSDP_FALSE;
113 DSDPFunctionReturn(0);
114}
115#undef __FUNCT__
116#define __FUNCT__ "DSDPSetDefaultParameters"
123
124 int info;
125 DSDPFunctionBegin;
126 DSDPValid(dsdp);
127
128 /* Stopping parameters */
129 info=DSDPSetMaxIts(dsdp,500);DSDPCHKERR(info);
130 info=DSDPSetGapTolerance(dsdp,1.0e-6);DSDPCHKERR(info);
131 info=DSDPSetPNormTolerance(dsdp,1.0e30);DSDPCHKERR(info);
132 if (dsdp->m<100){info=DSDPSetGapTolerance(dsdp,1.0e-7);DSDPCHKERR(info);}
133 if (dsdp->m>3000){info=DSDPSetGapTolerance(dsdp,5.0e-6);DSDPCHKERR(info);}
134 info=RConeSetType(dsdp->rcone,DSDPInfeasible);DSDPCHKERR(info);
135 info=DSDPSetDualBound(dsdp,1.0e20);DSDPCHKERR(info);
136 info=DSDPSetStepTolerance(dsdp,5.0e-2);DSDPCHKERR(info);
137 info=DSDPSetRTolerance(dsdp,1.0e-6);DSDPCHKERR(info);
138 info=DSDPSetPTolerance(dsdp,1.0e-4);DSDPCHKERR(info);
139 /* Solver options */
140 info=DSDPSetMaxTrustRadius(dsdp,1.0e10);DSDPCHKERR(info);
141 info=DSDPUsePenalty(dsdp,0);DSDPCHKERR(info);
142 info=DSDPSetInitialBarrierParameter(dsdp,-1.0);DSDPCHKERR(info);
143 info=DSDPSetPotentialParameter(dsdp,3.0);DSDPCHKERR(info);
144 info=DSDPUseDynamicRho(dsdp,1);DSDPCHKERR(info);
145 info=DSDPSetR0(dsdp,-1.0);DSDPCHKERR(info);
146 info=DSDPSetPenaltyParameter(dsdp,1.0e8);DSDPCHKERR(info);
147 info=DSDPReuseMatrix(dsdp,4);DSDPCHKERR(info);
148 if (dsdp->m>100){info=DSDPReuseMatrix(dsdp,7);DSDPCHKERR(info);}
149 if (dsdp->m>1000){info=DSDPReuseMatrix(dsdp,10);DSDPCHKERR(info);}
150 if (dsdp->m<=100){info=DSDPSetPotentialParameter(dsdp,5.0);DSDPCHKERR(info);DSDPCHKERR(info);}
151 dsdp->maxschurshift=1.0e-11;
152 dsdp->mu0=-1.0;
153 dsdp->slestype=2;
154 info = DSDPSetYBounds(dsdp,-1e7,1e7);DSDPCHKERR(info);
155 DSDPFunctionReturn(0);
156}
157
158#undef __FUNCT__
159#define __FUNCT__ "DSDPSetDefaultMonitors"
166
167 int info;
168
169 DSDPFunctionBegin;
170 DSDPValid(dsdp);
171 dsdp->nmonitors=0;
172 info=DSDPSetMonitor(dsdp,DSDPDefaultConvergence,(void*)&dsdp->conv); DSDPCHKERR(info);
173 DSDPFunctionReturn(0);
174}
175
191#undef __FUNCT__
192#define __FUNCT__ "DSDPSetUp"
193int DSDPSetup(DSDP dsdp){
194
195 int i,info;
196 DSDPFunctionBegin;
197 DSDPValid(dsdp);
198
199 /* Create the Work Vectors */
200 info = DSDPVecDuplicate(dsdp->y,&dsdp->rhs1);DSDPCHKERR(info);
201 info = DSDPVecDuplicate(dsdp->y,&dsdp->rhs2);DSDPCHKERR(info);
202 info = DSDPVecDuplicate(dsdp->y,&dsdp->rhs);DSDPCHKERR(info);
203 info = DSDPVecDuplicate(dsdp->y,&dsdp->rhstemp);DSDPCHKERR(info);
204 info = DSDPVecDuplicate(dsdp->y,&dsdp->dy1);DSDPCHKERR(info);
205 info = DSDPVecDuplicate(dsdp->y,&dsdp->dy2);DSDPCHKERR(info);
206 info = DSDPVecDuplicate(dsdp->y,&dsdp->dy);DSDPCHKERR(info);
207 info = DSDPVecDuplicate(dsdp->y,&dsdp->y0);DSDPCHKERR(info);
208 info = DSDPVecDuplicate(dsdp->y,&dsdp->xmakerrhs);DSDPCHKERR(info);
209 for (i=0;i<MAX_XMAKERS;i++){
210 info = DSDPVecDuplicate(dsdp->y,&dsdp->xmaker[i].y);DSDPCHKERR(info);
211 info = DSDPVecDuplicate(dsdp->y,&dsdp->xmaker[i].dy);DSDPCHKERR(info);
212 info = DSDPVecDuplicate(dsdp->y,&dsdp->xmaker[i].rhs);DSDPCHKERR(info);
213 }
214
215 /* Create M */
216 info = DSDPSetUpCones(dsdp);DSDPCHKERR(info);
217 info = DSDPSchurMatSetup(dsdp->M,dsdp->ytemp);DSDPCHKERR(info);
218
219 info = DSDPCGSetup(dsdp->sles,dsdp->ytemp); DSDPCHKERR(info);
220
221 info = DSDPSetUpCones2(dsdp,dsdp->y,dsdp->M);DSDPCHKERR(info);
222 info = DSDPGetConicDimension(dsdp,&dsdp->np);DSDPCHKERR(info);
223
224 info=DSDPComputeDataNorms(dsdp);DSDPCHKERR(info);
225 dsdp->pinfeas=dsdp->bnorm+1;
226 dsdp->perror=dsdp->bnorm+1;
227 info=DSDPScaleData(dsdp);DSDPCHKERR(info);
228
229 info=DSDPGetConicDimension(dsdp,&dsdp->np);DSDPCHKERR(info);
230 dsdp->solvetime=0;
231 dsdp->cgtime=0;
232 dsdp->ptime=0;
233 dsdp->dtime=0;
234 dsdp->ctime=0;
235 info=DSDPEventLogRegister("Primal Step",&dsdp->ptime);
236 info=DSDPEventLogRegister("Dual Step",&dsdp->dtime);
237 info=DSDPEventLogRegister("Corrector Step",&dsdp->ctime);
238 info=DSDPEventLogRegister("CG Solve",&dsdp->cgtime);
239 info=DSDPEventLogRegister("DSDP Solve",&dsdp->solvetime);
240 dsdp->setupcalled=DSDP_TRUE;
241 DSDPFunctionReturn(0);
242}
243
244
245
246#undef __FUNCT__
247#define __FUNCT__ "DSDPGetSchurMatrix"
248int DSDPGetSchurMatrix(DSDP dsdp, DSDPSchurMat *M){
249 DSDPFunctionBegin;
250 DSDPValid(dsdp);
251 *M=dsdp->M;
252 DSDPFunctionReturn(0);
253}
254
255#undef __FUNCT__
256#define __FUNCT__ "DSDPGetConvergenceMonitor"
268int DSDPGetConvergenceMonitor(DSDP dsdp, ConvergenceMonitor**ctx){
269 DSDPFunctionBegin;
270 DSDPValid(dsdp);
271 *ctx=&dsdp->conv;
272 DSDPFunctionReturn(0);
273}
274
275
276#undef __FUNCT__
277#define __FUNCT__ "DSDPComputeDataNorms"
284 int info;
285 DSDPVec ytemp=dsdp->ytemp;
286 DSDPFunctionBegin;
287 DSDPValid(dsdp);
288 info = DSDPComputeANorm2(dsdp,ytemp);DSDPCHKERR(info);
289 info = DSDPFixedVariablesNorm(dsdp->M,ytemp);DSDPCHKERR(info);
290 info = DSDPVecGetC(ytemp,&dsdp->cnorm);DSDPCHKERR(info);
291 dsdp->cnorm=sqrt(dsdp->cnorm);
292 info = DSDPVecSetR(ytemp,0);DSDPCHKERR(info);
293 info = DSDPVecSetC(ytemp,0);DSDPCHKERR(info);
294 info = DSDPVecNorm1(ytemp,&dsdp->anorm);DSDPCHKERR(info);
295 dsdp->anorm=sqrt(dsdp->anorm);
296 DSDPLogInfo(0,2,"Norm of data: %4.2e\n",dsdp->anorm);
297 info=DSDPVecCopy(dsdp->b,ytemp);DSDPCHKERR(info);
298 info = DSDPVecSetR(ytemp,0);DSDPCHKERR(info);
299 info = DSDPVecSetC(ytemp,0);DSDPCHKERR(info);
300 info = DSDPVecNorm2(ytemp,&dsdp->bnorm);DSDPCHKERR(info);
301 DSDPFunctionReturn(0);
302}
303
304#undef __FUNCT__
305#define __FUNCT__ "DSDPScaleData"
312 int info;
313 double scale;
314 DSDPFunctionBegin;
315 DSDPValid(dsdp);
316 scale=1.0*dsdp->anorm;
317 if (dsdp->bnorm){ scale/=dsdp->bnorm;}
318 if (dsdp->cnorm){ scale/=dsdp->cnorm;}
319 scale=DSDPMin(scale,1.0);
320 scale=DSDPMax(scale,1.0e-6);
321 if (dsdp->cnorm==0){ scale=1;}
322 info=DSDPSetScale(dsdp,scale);DSDPCHKERR(info);
323 DSDPFunctionReturn(0);
324}
325
341#undef __FUNCT__
342#define __FUNCT__ "DSDPSolve"
343int DSDPSolve(DSDP dsdp){
344 int info;
345 DSDPFunctionBegin;
346 info=DSDPEventLogBegin(dsdp->solvetime);
347 dsdp->pdfeasible=DSDP_PDUNKNOWN;
348 info=DSDPSetConvergenceFlag(dsdp,CONTINUE_ITERATING);DSDPCHKERR(info);
349 info=DSDPInitializeVariables(dsdp);DSDPCHKERR(info);
350 info=DSDPSolveDynamicRho(dsdp);DSDPCHKERR(info);
351 if (dsdp->pstep==1){info=DSDPRefineStepDirection(dsdp,dsdp->xmakerrhs,dsdp->xmaker[0].dy);DSDPCHKERR(info);}
352 if (dsdp->pdfeasible==DSDP_PDUNKNOWN) dsdp->pdfeasible=DSDP_PDFEASIBLE;
353 info=DSDPEventLogEnd(dsdp->solvetime);
354 DSDPFunctionReturn(0);
355}
356
357
358#undef __FUNCT__
359#define __FUNCT__ "DSDPCallMonitors"
367int DSDPCallMonitors(DSDP dsdp,DMonitor dmonitor[], int ndmonitors){
368 int i,info;
369 DSDPFunctionBegin;
370 for (i=0; i<ndmonitors;i++){
371 info=(dmonitor[i].monitor)(dsdp,dmonitor[i].monitorctx); DSDPCHKERR(info);
372 }
373 DSDPFunctionReturn(0);
374}
375/* ---------------------------------------------------------- */
376#undef __FUNCT__
377#define __FUNCT__ "DSDPCheckConvergence"
385 int info;
386 DSDPTruth unbounded;
387
388 DSDPFunctionBegin;
389 info = DSDPGetConicDimension(dsdp,&dsdp->np);DSDPCHKERR(info);
390 dsdp->rgap=(dsdp->ppobj-dsdp->ddobj)/(1.0+fabs(dsdp->ppobj)+fabs(dsdp->ddobj));
391 dsdp->pstepold=dsdp->pstep;
392 if (dsdp->reason==CONTINUE_ITERATING){
393 if (dsdp->itnow>0){
394 info=DSDPCheckForUnboundedObjective(dsdp,&unbounded);DSDPCHKERR(info);
395 if (unbounded==DSDP_TRUE){
396 dsdp->pdfeasible=DSDP_UNBOUNDED;
397 info=DSDPSetConvergenceFlag(dsdp,DSDP_CONVERGED); DSDPCHKERR(info);
398 }
399 }
400 if (dsdp->reason==CONTINUE_ITERATING){
401 if (dsdp->muold<dsdp->mutarget && dsdp->pstep==1 && dsdp->dstep==1 && dsdp->rgap<1e-5){
402 info=DSDPSetConvergenceFlag(dsdp,DSDP_NUMERICAL_ERROR); DSDPCHKERR(info);
403 DSDPLogInfo(0,2,"DSDP Finished: Numerical issues: Increase in Barrier function. \n");}
404 if (dsdp->itnow >= dsdp->maxiter){
405 info=DSDPSetConvergenceFlag(dsdp,DSDP_MAX_IT); DSDPCHKERR(info);}
406 if (dsdp->Mshift>dsdp->maxschurshift){
407 info = DSDPSetConvergenceFlag(dsdp,DSDP_INDEFINITE_SCHUR_MATRIX); DSDPCHKERR(info);
408 }
409 }
410 info=DSDPCallMonitors(dsdp,dsdp->dmonitor,dsdp->nmonitors);DSDPCHKERR(info);
411 info=DSDPMonitorCones(dsdp,0); DSDPCHKERR(info);
412 }
413 dsdp->muold=dsdp->mutarget;
414 info = DSDPStopReason(dsdp,reason); DSDPCHKERR(info);
415 DSDPFunctionReturn(0);
416}
417
418
419
420/* ---------------------------------------------------------- */
421#undef __FUNCT__
422#define __FUNCT__ "DSDPTakeDown"
429
430 int i,info;
431
432 DSDPFunctionBegin;
433 DSDPValid(dsdp);
434 info = DSDPVecDestroy(&dsdp->rhs);DSDPCHKERR(info);
435 info = DSDPVecDestroy(&dsdp->rhs1);DSDPCHKERR(info);
436 info = DSDPVecDestroy(&dsdp->rhs2);DSDPCHKERR(info);
437 info = DSDPVecDestroy(&dsdp->rhstemp);DSDPCHKERR(info);
438 info = DSDPVecDestroy(&dsdp->y);DSDPCHKERR(info);
439 info = DSDPVecDestroy(&dsdp->ytemp);DSDPCHKERR(info);
440 info = DSDPVecDestroy(&dsdp->dy1);DSDPCHKERR(info);
441 info = DSDPVecDestroy(&dsdp->dy2);DSDPCHKERR(info);
442 info = DSDPVecDestroy(&dsdp->dy);DSDPCHKERR(info);
443 for (i=0;i<MAX_XMAKERS;i++){
444 info = DSDPVecDestroy(&dsdp->xmaker[i].y);DSDPCHKERR(info);
445 info = DSDPVecDestroy(&dsdp->xmaker[i].dy);DSDPCHKERR(info);
446 info = DSDPVecDestroy(&dsdp->xmaker[i].rhs);DSDPCHKERR(info);
447 }
448 info = DSDPVecDestroy(&dsdp->xmakerrhs);DSDPCHKERR(info);
449 info = DSDPVecDestroy(&dsdp->y0);DSDPCHKERR(info);
450 info = DSDPVecDestroy(&dsdp->b);DSDPCHKERR(info);
451
452 info = DSDPCGDestroy(&dsdp->sles);DSDPCHKERR(info);
453 info = DSDPDestroyCones(dsdp);DSDPCHKERR(info);
454 info = DSDPSchurMatDestroy(&dsdp->M);DSDPCHKERR(info);
455 info = DSDPGetConicDimension(dsdp,&dsdp->np);DSDPCHKERR(info);
456 dsdp->setupcalled=DSDP_FALSE;
457 DSDPFunctionReturn(0);
458}
459
469int DSDPSetDestroyRoutine(DSDP dsdp, int (*fd)(void*), void* ctx){
470 int nd=dsdp->ndroutines;
471 if (nd<10){
472 dsdp->droutine[nd].f=fd;
473 dsdp->droutine[nd].ptr=ctx;
474 dsdp->ndroutines++;
475 } else {
476 printf("TOO MANY Destroy routines\n");
477 return 1;
478 }
479 return 0;
480}
481
482
494#undef __FUNCT__
495#define __FUNCT__ "DSDPDestroy"
497 int i,info;
498 DSDPFunctionBegin;
499 DSDPValid(dsdp);
500 for (i=0;i<dsdp->ndroutines;i++){
501 info=(*dsdp->droutine[i].f)(dsdp->droutine[i].ptr);DSDPCHKERR(info);
502 }
503 info=DSDPTakeDown(dsdp);DSDPCHKERR(info);
504 DSDPFREE(&dsdp,&info);DSDPCHKERR(info);
505 DSDPFunctionReturn(0);
506}
int DSDPCreateLUBoundsCone(DSDP dsdp, LUBounds *dspcone)
Create bounds cone.
Definition allbounds.c:566
The API to DSDP for those applications using DSDP as a subroutine library.
Internal data structure for the DSDP solver.
int DSDPSetUpCones(DSDP)
Each cone should factor data or allocate internal data structures.
Definition dsdpcops.c:58
int DSDPComputeANorm2(DSDP, DSDPVec)
Compute norm of A and C.
Definition dsdpcops.c:246
int DSDPMonitorCones(DSDP, int)
This routine is called once per iteration.
Definition dsdpcops.c:450
int DSDPSolveDynamicRho(DSDP)
Apply dual-scaling algorithm.
Definition dualalg.c:121
int DSDPGetConicDimension(DSDP, double *)
Get the total dimension of the cones.
Definition dsdpcops.c:401
int DSDPDestroyCones(DSDP)
Each cone shoudl free its data structures.
Definition dsdpcops.c:107
int DSDPInitializeVariables(DSDP)
Initialize variables and factor S.
Definition dualalg.c:475
int DSDPSetUpCones2(DSDP, DSDPVec, DSDPSchurMat)
Each cone should allocate its data structures .
Definition dsdpcops.c:84
DSDPTruth
Boolean variables.
@ DSDP_FALSE
@ DSDP_TRUE
struct DSDP_C * DSDP
An implementation of the dual-scaling algorithm for semidefinite programming.
DSDPTerminationReason
There are many reasons to terminate the solver.
@ DSDP_INDEFINITE_SCHUR_MATRIX
@ CONTINUE_ITERATING
@ DSDP_MAX_IT
@ DSDP_CONVERGED
@ DSDP_NUMERICAL_ERROR
@ DSDP_UNBOUNDED
@ DSDP_PDFEASIBLE
@ DSDP_PDUNKNOWN
int DSDPSchurMatSetup(DSDPSchurMat M, DSDPVec Y)
Set up the data structure.
int DSDPSchurMatInitialize(DSDPSchurMat *M)
Initialize pointers to null.
int DSDPSchurMatDestroy(DSDPSchurMat *M)
Free the memory in the data structure.
struct DSDPSchurMat_C DSDPSchurMat
This object represents the Schur Matrix. Its structure is opaque to the DSDP solver,...
int DSDPSetDefaultParameters(DSDP dsdp)
Set default parameters.
Definition dsdpsetup.c:122
int DSDPScaleData(DSDP dsdp)
Scale the matrix C.
Definition dsdpsetup.c:311
int DSDPTakeDown(DSDP dsdp)
Destroy internal data structures.
Definition dsdpsetup.c:428
int DSDPCheckConvergence(DSDP dsdp, DSDPTerminationReason *reason)
Check for convergence and monitor solution.
Definition dsdpsetup.c:384
int DSDPCallMonitors(DSDP dsdp, DMonitor dmonitor[], int ndmonitors)
Call the monitor routines.
Definition dsdpsetup.c:367
int DSDPGetConvergenceMonitor(DSDP dsdp, ConvergenceMonitor **ctx)
Get the structure containing convergence parameters.
Definition dsdpsetup.c:268
int DSDPSetDefaultMonitors(DSDP dsdp)
Set convergence monitor.
Definition dsdpsetup.c:165
int DSDPComputeDataNorms(DSDP dsdp)
Compute norms of A,C, and b.
Definition dsdpsetup.c:283
int DSDPSetDefaultStatistics(DSDP dsdp)
Set default statistics.
Definition dsdpsetup.c:84
int DSDPSetDestroyRoutine(DSDP dsdp, int(*fd)(void *), void *ctx)
Set a routine that will be called during DSDPDestroy().
Definition dsdpsetup.c:469
Error handling, printing, and profiling.
struct DSDPVec_C DSDPVec
This object hold m+2 variables: a scaling of C, the y variables, and r.
Definition dsdpvec.h:25
int DSDPSetup(DSDP dsdp)
Set up data structures in the solver and the cones associated with it.
Definition dsdpsetup.c:193
int DSDPDestroy(DSDP dsdp)
Free the internal data structures of the solver and the cones associated with it.
Definition dsdpsetup.c:496
int DSDPCreate(int m, DSDP *dsdpnew)
Create a DSDP solver. FIRST DSDP routine!
Definition dsdpsetup.c:30
int DSDPSolve(DSDP dsdp)
Apply DSDP to the problem.
Definition dsdpsetup.c:343
int DSDPDefaultConvergence(DSDP, void *)
Check for Convergence.
int DSDPSetConvergenceFlag(DSDP dsdp, DSDPTerminationReason reason)
Monitor each iteration of the solver.
int DSDPSetGapTolerance(DSDP dsdp, double gaptol)
Terminate the solver when the relative duality gap is less than this tolerance.
int DSDPSetPNormTolerance(DSDP dsdp, double ptol)
Terminate the solver when the relative duality gap is suffiently small and the PNorm is less than thi...
int DSDPStopReason(DSDP dsdp, DSDPTerminationReason *reason)
Copy the reason why the solver terminated.
int DSDPSetStepTolerance(DSDP dsdp, double steptol)
Terminate the solver if the step length in (DD) is below this tolerance.
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 DSDPSetMaxIts(DSDP dsdp, int its)
Terminate the solver after this number of iterations.
int DSDPSetDualBound(DSDP dsdp, double dbound)
Terminate the solver if the objective value in (DD) is greater than this tolerance.
int DSDPSetPTolerance(DSDP dsdp, double inftol)
Classify (P) as feasible only if the infeasibility is less than this tolerance.
Definition dsdpx.c:365
int DSDPSetR0(DSDP dsdp, double res)
Set an initial value for the variable r in (DD)
int DSDPUsePenalty(DSDP dsdp, int yesorno)
Use penalty parameter to enforce feasibility.
int DSDPSetPotentialParameter(DSDP dsdp, double rho)
Set the potential parameter.
int DSDPSetMonitor(DSDP dsdp, int(*monitor)(DSDP, void *), void *monitorctx)
Monitor each iteration of the solver.
int DSDPReuseMatrix(DSDP dsdp, int rm)
Reuse the Hessian of the barrier function multiple times at each DSDP iteration.
int DSDPSetPenaltyParameter(DSDP dsdp, double Gamma)
Set the penalty parameter Gamma.
int DSDPSetScale(DSDP dsdp, double scale)
Set the internal scaling factor.
int DSDPSetMaxTrustRadius(DSDP dsdp, double rad)
Set a maximum trust radius on the step direction.
int DSDPSetYBounds(DSDP dsdp, double lbound, double ubound)
Bound the variables y.
int DSDPUseDynamicRho(DSDP dsdp, int yesorno)
Use a dynamic strategy to choose parameter rho.