10 #define __FUNCT__ "DSDPInspectXY"
11 int DSDPInspectXY(
DSDP dsdp,
double xmakermu,
DSDPVec xmakery,
DSDPVec xmakerdy,
DSDPVec AX,
double *tracexs2,
double *pobj2,
double *rpinfeas2){
15 info=BoundYConeAddX(dsdp->ybcone,xmakermu,xmakery,xmakerdy,AX,tracexs2); DSDPCHKERR(info);
16 info=DSDPVecGetC(AX,pobj2);DSDPCHKERR(info);
18 info=DSDPVecSetC(AX,0);DSDPCHKERR(info);
19 info=DSDPVecSetR(AX,0);DSDPCHKERR(info);
20 info=DSDPVecNorm1(AX,rpinfeas2); DSDPCHKERR(info);
21 DSDPFunctionReturn(0);
54 #define __FUNCT__ "DSDPComputeX"
57 double pobj=0,ppobj2=0,ddobj,tracexs=0,tracexs2=0,rpinfeas=0,rpinfeas2=0,rpobjerr=0;
58 double err1,cc,rrr,bigM,ymax,pfeastol=dsdp->pinfeastol;
66 info=
DSDPGetR(dsdp,&rrr); DSDPCHKERR(info);
67 info=DSDPGetPenalty(dsdp,&bigM);DSDPCHKERR(info);
71 for (i=0;i<MAX_XMAKERS;i++){
72 if (i>0 && dsdp->xmaker[i].pstep<1)
continue;
73 info=
DSDPComputeXVariables(dsdp,dsdp->xmaker[i].mu,dsdp->xmaker[i].y,dsdp->xmaker[i].dy,AX,&tracexs);DSDPCHKERR(info);
74 info=DSDPVecGetC(AX,&pobj); DSDPCHKERR(info);
75 info=DSDPVecGetR(AX,&dsdp->tracex); DSDPCHKERR(info);
76 info=DSDPVecSetC(AX,0);DSDPCHKERR(info);
77 info=DSDPVecSetR(AX,0);DSDPCHKERR(info);
78 info=DSDPVecNormInfinity(AX,&rpinfeas);DSDPCHKERR(info);
79 rpinfeas=rpinfeas/(dsdp->tracex+1);
81 DSDPLogInfo(0,2,
"POBJ: %4.4e, DOBJ: %4.4e\n",pobj,ddobj/cc);
83 info=DSDPVecNorm2(AX,&err1);DSDPCHKERR(info);
84 dsdp->tracexs=tracexs;
88 info=DSDPInspectXY(dsdp,dsdp->xmaker[i].mu,dsdp->xmaker[i].y,dsdp->xmaker[i].dy,AX,&tracexs2,&ppobj2,&rpinfeas2);DSDPCHKERR(info);
89 rpinfeas2=rpinfeas2/(dsdp->tracex+1);
92 DSDPLogInfo(0,2,
"X P Infeas: %4.2e , PObj: %4.8e\n",rpinfeas,pobj*(cc));
93 DSDPLogInfo(0,2,
"TOTAL P Infeas: %4.2e PObj: %4.8e\n",rpinfeas2,ppobj2*(cc));
94 rpobjerr= fabs(pobj-dsdp->ppobj)/(1+fabs(dsdp->ppobj));
96 if (rpinfeas2 < pfeastol){
99 if (rpinfeas>pfeastol/100 && fabs(rrr)>dsdp->dinfeastol){
101 DSDPLogInfo(0,2,
"Warning: Try Increasing penalty parameter\n");
102 }
else if (rpinfeas>pfeastol && ddobj>0 && ppobj2<0 && fabs(rrr)<dsdp->dinfeastol){
104 DSDPLogInfo(0,2,
"Warning: D probably unbounded\n");
106 }
else if ( fabs(rrr)>dsdp->dinfeastol){
108 DSDPLogInfo(0,2,
"Warning: D probably infeasible \n");
116 DSDPLogInfo(0,2,
"Try backup X\n");
122 DSDPFunctionReturn(0);
127 #define __FUNCT__ "DSDPSaveBackupYForX"
128 int DSDPSaveBackupYForX(
DSDP dsdp,
int count,
double mu,
double pstep){
132 info=DSDPVecCopy(dsdp->y,dsdp->xmaker[count].y);DSDPCHKERR(info);
133 info=
DSDPComputeDY(dsdp,mu,dsdp->xmaker[count].dy,&pnorm); DSDPCHKERR(info);
134 dsdp->xmaker[count].pstep=pstep;
135 dsdp->xmaker[count].mu = mu;
136 DSDPFunctionReturn(0);
148 #define __FUNCT__ "DSDPSaveYForX"
151 double pnorm,newgap,ymax,dd=0;
153 dsdp->pstepold=dsdp->pstep;
156 info=DSDPVecCopy(dsdp->y,dsdp->xmaker[0].y);DSDPCHKERR(info);
157 dsdp->xmaker[0].pstep=pstep;
158 }
else if (dsdp->Mshift*ymax>dsdp->pinfeastol*10){
160 if (pstep==1 && newgap>0){
161 dsdp->ppobj = dsdp->ddobj + newgap; dsdp->mu=(newgap)/(dsdp->np);
162 dsdp->dualitygap=newgap;
164 info=DSDPVecZero(dsdp->rhstemp); DSDPCHKERR(info);
165 info=BoundYConeAddX(dsdp->ybcone,dsdp->xmaker[0].mu,dsdp->xmaker[0].y,dsdp->xmaker[0].dy,dsdp->rhstemp,&dd); DSDPCHKERR(info);
166 info=DSDPVecSetC(dsdp->rhstemp,0);
167 info=DSDPVecSetR(dsdp->rhstemp,0);
168 info=DSDPVecNormInfinity(dsdp->rhstemp,&dsdp->pinfeas); DSDPCHKERR(info);
169 dsdp->pinfeas+=dsdp->Mshift*ymax;
170 if (0==1){info=DSDPVecView(dsdp->rhstemp);}
173 info=DSDPVecCopy(dsdp->y,dsdp->xmaker[0].y);DSDPCHKERR(info);
174 dsdp->xmaker[0].pstep=pstep;
176 info=
DSDPComputeDY(dsdp,mu,dsdp->xmaker[0].dy,&pnorm); DSDPCHKERR(info);
177 dsdp->xmaker[0].mu = mu;
179 if (pstep==1 && newgap>0){
180 dsdp->ppobj = dsdp->ddobj + newgap; dsdp->mu=(newgap)/(dsdp->np);
181 dsdp->dualitygap=newgap;
183 info=DSDPVecZero(dsdp->rhstemp); DSDPCHKERR(info);
184 info=BoundYConeAddX(dsdp->ybcone,dsdp->xmaker[0].mu,dsdp->xmaker[0].y,dsdp->xmaker[0].dy,dsdp->rhstemp,&dd); DSDPCHKERR(info);
185 info=DSDPVecSetC(dsdp->rhstemp,0);
186 info=DSDPVecSetR(dsdp->rhstemp,0);
187 info=DSDPVecNormInfinity(dsdp->rhstemp,&dsdp->pinfeas); DSDPCHKERR(info);
188 dsdp->pinfeas+=dsdp->Mshift*ymax;
189 if (0==1){info=DSDPVecView(dsdp->rhstemp);}
197 info=
DSDPPassXVectors(dsdp,dsdp->xmaker[0].mu,dsdp->xmaker[0].y,dsdp->xmaker[0].dy); DSDPCHKERR(info);
198 info=
DSDPGetRR(dsdp,&r);DSDPCHKERR(info);
199 if (r&& dsdp->rgap<0.1){
200 info=RConeGetRX(dsdp->rcone,&rx);DSDPCHKERR(info);
201 info=DSDPGetPenalty(dsdp,&penalty);DSDPCHKERR(info);
202 dsdp->pinfeas = dsdp->pinfeas *(1+fabs(penalty-rx));
206 if (pstep==1.0 && dsdp->rgap>5.0e-1) {
207 info=DSDPSaveBackupYForX(dsdp,MAX_XMAKERS-1,mu,pstep);DSDPCHKERR(info);
209 if (pstep==1.0 && dsdp->rgap>1.0e-3) {
210 info=DSDPSaveBackupYForX(dsdp,2,mu,pstep);DSDPCHKERR(info);
212 if (pstep==1.0 && dsdp->rgap>1.0e-5) {
213 info=DSDPSaveBackupYForX(dsdp,1,mu,pstep);DSDPCHKERR(info);
216 DSDPFunctionReturn(0);
231 #define __FUNCT__ "DSDPGetPObjective"
238 *pobj=(dsdp->pobj)/scale;
239 DSDPFunctionReturn(0);
253 #define __FUNCT__ "DSDPGetSolutionType"
257 *pdfeasible=dsdp->pdfeasible;
258 DSDPFunctionReturn(0);
277 #define __FUNCT__ "DSDPGetTraceX"
281 *tracex=dsdp->tracex;
282 DSDPFunctionReturn(0);
296 #define __FUNCT__ "DSDPGetFinalErrors"
299 double scale,rr,bnorm,dobj=0,pobj=0;
303 info=DSDPVecGetR(dsdp->y,&rr); DSDPCHKERR(info);
308 err[2]=fabs(rr)/scale;
311 err[5]=dsdp->tracexs/scale;
313 err[2] /= (1.0+dsdp->cnorm);
314 info=DSDPVecCopy(dsdp->b,dsdp->ytemp);DSDPCHKERR(info);
315 info=DSDPVecSetC(dsdp->ytemp,0);DSDPCHKERR(info);
316 info=DSDPVecSetR(dsdp->ytemp,0);DSDPCHKERR(info);
317 info=DSDPVecNormInfinity(dsdp->ytemp,&bnorm);
318 err[0]=dsdp->perror/(1.0+bnorm);
320 err[4]=(err[4])/(1.0+fabs(pobj)+fabs(dobj));
321 err[5]=(err[5])/(1.0+fabs(pobj)+fabs(dobj));
322 DSDPFunctionReturn(0);
342 #define __FUNCT__ "DSDPGetPInfeasibility"
346 if (pperror) *pperror=dsdp->pinfeas;
347 DSDPFunctionReturn(0);
364 #define __FUNCT__ "DSDPSetPTolerance"
368 if (inftol > 0) dsdp->pinfeastol = inftol;
369 DSDPLogInfo(0,2,
"Set P Infeasibility Tolerance: %4.4e\n",inftol);
370 DSDPFunctionReturn(0);
385 #define __FUNCT__ "DSDPGetPTolerance"
389 if (inftol) *inftol=dsdp->pinfeastol;
390 DSDPFunctionReturn(0);
408 #define __FUNCT__ "DSDPSetRTolerance"
412 if (inftol > 0) dsdp->dinfeastol = inftol;
413 DSDPLogInfo(0,2,
"Set D Infeasibility Tolerance: %4.4e\n",inftol);
414 DSDPFunctionReturn(0);
433 #define __FUNCT__ "DSDPGetRTolerance"
437 *inftol=dsdp->dinfeastol;
438 DSDPFunctionReturn(0);
454 #define __FUNCT__ "DSDPGetYMakeX"
460 if (dsdp->m < m-1) DSDPFunctionReturn(1);
461 if (dsdp->m > m) DSDPFunctionReturn(1);
462 info=DSDPVecCopy(dsdp->xmaker[0].y,dsdp->ytemp); DSDPCHKERR(info);
464 info=DSDPVecGetArray(dsdp->ytemp,&yy);DSDPCHKERR(info);
465 for (i=0;i<m;i++) y[i]=yy[i+1]/scale;
466 info=DSDPVecRestoreArray(dsdp->ytemp,&yy);DSDPCHKERR(info);
467 DSDPFunctionReturn(0);
482 #define __FUNCT__ "DSDPGetDYMakeX"
488 if (dsdp->m < m-1) DSDPFunctionReturn(1);
489 if (dsdp->m > m) DSDPFunctionReturn(1);
490 info=DSDPVecCopy(dsdp->xmaker[0].dy,dsdp->ytemp); DSDPCHKERR(info);
492 info=DSDPVecGetArray(dsdp->ytemp,&yy);DSDPCHKERR(info);
493 for (i=0;i<m;i++) dy[i]=yy[i+1]/scale;
494 info=DSDPVecRestoreArray(dsdp->ytemp,&yy);DSDPCHKERR(info);
495 DSDPFunctionReturn(0);
510 #define __FUNCT__ "DSDPGetMuMakeX"
517 *mu=dsdp->xmaker[0].mu/scale;
518 DSDPFunctionReturn(0);