29#define __FUNCT__ "DSDPCreate"
30int DSDPCreate(
int m,
DSDP* dsdpnew){
37 DSDPCALLOC1(&dsdp,
PD_DSDP,&info);DSDPCHKERR(info);
42 DSDPEventLogInitialize();
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);
58 info = DSDPAddRCone(dsdp,&dsdp->rcone);DSDPCHKERR(info);
67 info = DSDPSetDefaultSchurMatrixStructure(dsdp); DSDPCHKERR(info);
68 info = DSDPCGInitialize(&dsdp->sles); DSDPCHKERR(info);
73 DSDPFunctionReturn(0);
78#define __FUNCT__ "DSDPSetDefaultStatistics"
96 dsdp->dualitygap=dsdp->ppobj-dsdp->ddobj;
99 for (i=0;i<MAX_XMAKERS;i++){
100 dsdp->xmaker[i].mu=1.0e200;
101 dsdp->xmaker[i].pstep=0.0;
113 DSDPFunctionReturn(0);
116#define __FUNCT__ "DSDPSetDefaultParameters"
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);
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;
154 info = DSDPSetYBounds(dsdp,-1e7,1e7);DSDPCHKERR(info);
155 DSDPFunctionReturn(0);
159#define __FUNCT__ "DSDPSetDefaultMonitors"
172 info=DSDPSetMonitor(dsdp,DSDPDefaultConvergence,(
void*)&dsdp->conv); DSDPCHKERR(info);
173 DSDPFunctionReturn(0);
192#define __FUNCT__ "DSDPSetUp"
193int DSDPSetup(
DSDP dsdp){
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);
219 info = DSDPCGSetup(dsdp->sles,dsdp->ytemp); DSDPCHKERR(info);
225 dsdp->pinfeas=dsdp->bnorm+1;
226 dsdp->perror=dsdp->bnorm+1;
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);
241 DSDPFunctionReturn(0);
247#define __FUNCT__ "DSDPGetSchurMatrix"
252 DSDPFunctionReturn(0);
256#define __FUNCT__ "DSDPGetConvergenceMonitor"
268int DSDPGetConvergenceMonitor(
DSDP dsdp, ConvergenceMonitor**ctx){
272 DSDPFunctionReturn(0);
277#define __FUNCT__ "DSDPComputeDataNorms"
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);
305#define __FUNCT__ "DSDPScaleData"
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);
342#define __FUNCT__ "DSDPSolve"
343int DSDPSolve(
DSDP dsdp){
346 info=DSDPEventLogBegin(dsdp->solvetime);
351 if (dsdp->pstep==1){info=DSDPRefineStepDirection(dsdp,dsdp->xmakerrhs,dsdp->xmaker[0].dy);DSDPCHKERR(info);}
353 info=DSDPEventLogEnd(dsdp->solvetime);
354 DSDPFunctionReturn(0);
359#define __FUNCT__ "DSDPCallMonitors"
370 for (i=0; i<ndmonitors;i++){
371 info=(dmonitor[i].monitor)(dsdp,dmonitor[i].monitorctx); DSDPCHKERR(info);
373 DSDPFunctionReturn(0);
377#define __FUNCT__ "DSDPCheckConvergence"
390 dsdp->rgap=(dsdp->ppobj-dsdp->ddobj)/(1.0+fabs(dsdp->ppobj)+fabs(dsdp->ddobj));
391 dsdp->pstepold=dsdp->pstep;
394 info=DSDPCheckForUnboundedObjective(dsdp,&unbounded);DSDPCHKERR(info);
397 info=DSDPSetConvergenceFlag(dsdp,
DSDP_CONVERGED); DSDPCHKERR(info);
401 if (dsdp->muold<dsdp->mutarget && dsdp->pstep==1 && dsdp->dstep==1 && dsdp->rgap<1e-5){
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){
410 info=
DSDPCallMonitors(dsdp,dsdp->dmonitor,dsdp->nmonitors);DSDPCHKERR(info);
413 dsdp->muold=dsdp->mutarget;
414 info = DSDPStopReason(dsdp,reason); DSDPCHKERR(info);
415 DSDPFunctionReturn(0);
422#define __FUNCT__ "DSDPTakeDown"
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);
448 info = DSDPVecDestroy(&dsdp->xmakerrhs);DSDPCHKERR(info);
449 info = DSDPVecDestroy(&dsdp->y0);DSDPCHKERR(info);
450 info = DSDPVecDestroy(&dsdp->b);DSDPCHKERR(info);
452 info = DSDPCGDestroy(&dsdp->sles);DSDPCHKERR(info);
457 DSDPFunctionReturn(0);
470 int nd=dsdp->ndroutines;
472 dsdp->droutine[nd].f=fd;
473 dsdp->droutine[nd].ptr=ctx;
476 printf(
"TOO MANY Destroy routines\n");
495#define __FUNCT__ "DSDPDestroy"
496int DSDPDestroy(
DSDP dsdp){
500 for (i=0;i<dsdp->ndroutines;i++){
501 info=(*dsdp->droutine[i].f)(dsdp->droutine[i].ptr);DSDPCHKERR(info);
504 DSDPFREE(&dsdp,&info);DSDPCHKERR(info);
505 DSDPFunctionReturn(0);
int DSDPCreateLUBoundsCone(DSDP dsdp, LUBounds *dspcone)
Create bounds cone.
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.
int DSDPComputeANorm2(DSDP, DSDPVec)
Compute norm of A and C.
int DSDPMonitorCones(DSDP, int)
This routine is called once per iteration.
int DSDPSolveDynamicRho(DSDP)
Apply dual-scaling algorithm.
int DSDPGetConicDimension(DSDP, double *)
Get the total dimension of the cones.
int DSDPDestroyCones(DSDP)
Each cone shoudl free its data structures.
int DSDPInitializeVariables(DSDP)
Initialize variables and factor S.
int DSDPSetUpCones2(DSDP, DSDPVec, DSDPSchurMat)
Each cone should allocate its data structures .
DSDPTerminationReason
There are many reasons to terminate the solver.
@ DSDP_INDEFINITE_SCHUR_MATRIX
DSDPTruth
Boolean variables.
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.
int DSDPSetDefaultParameters(DSDP dsdp)
Set default parameters.
int DSDPScaleData(DSDP dsdp)
Scale the matrix C.
int DSDPTakeDown(DSDP dsdp)
Destroy internal data structures.
int DSDPCheckConvergence(DSDP dsdp, DSDPTerminationReason *reason)
Check for convergence and monitor solution.
int DSDPCallMonitors(DSDP dsdp, DMonitor dmonitor[], int ndmonitors)
Call the monitor routines.
int DSDPSetDefaultMonitors(DSDP dsdp)
Set convergence monitor.
int DSDPComputeDataNorms(DSDP dsdp)
Compute norms of A,C, and b.
int DSDPSetDefaultStatistics(DSDP dsdp)
Set default statistics.
int DSDPSetDestroyRoutine(DSDP dsdp, int(*fd)(void *), void *ctx)
Set a routine that will be called during DSDPDestroy().
Error handling, printing, and profiling.
struct DSDPVec_C DSDPVec
This object hold m+2 variables: a scaling of C, the y variables, and r.
Schur complement matrix whose solution is the Newton direction.
Internal structures for the DSDP solver.