DSDP
dsdpstep.c
Go to the documentation of this file.
1#include "dsdpdualmat.h"
2#include "dsdpdsmat.h"
3#include "dsdpxmat.h"
4#include "dsdpsys.h"
5#include "dsdplanczos.h"
6#include "dsdplapack.h"
7
13typedef struct _P_Mat3* Mat3;
14
15static int MatMult3(Mat3 A, SDPConeVec x, SDPConeVec y);
16static int ComputeStepROBUST(Mat3 A, SDPConeVec *Q, int m, SDPConeVec W, SDPConeVec R,double*, SDPConeVec QAQTv, double *dwork, double *maxstep, double *mineig);
17static int ComputeStepFAST(Mat3 A, SDPConeVec *Q, int m, SDPConeVec W, double *dwork, int*iwork,double *maxstep, double *mineig);
18
19extern int DSDPGetEigsSTEP(double[],int,double[],int,long int[],int,
20 double[],int,double[],int,int[],int);
21
22int DSDPGetTriDiagonalEigs(int n,double *D,double *E, double*WORK2N,int*);
23
24struct _P_Mat3{
25 int type;
26 DSDPDualMat ss;
27 DSDPDSMat ds;
28 SDPConeVec V;
29 DSDPVMat x;
30};
31
32
33int DSDPGetTriDiagonalEigs(int N,double D[],double E[], double WORK[],int IIWORK[]){
34 ffinteger LDZ=DSDPMax(1,N),INFO,NN=N;
35 ffinteger M,IL=N-1,IU=N,*ISUPPZ=0;
36 ffinteger *IWORK=(ffinteger*)IIWORK,LIWORK,LWORK;
37 double WW[2],VL=-1e10,VU=1e10,*Z=0,ABSTOL=0;
38 char JOBZ='N', RANGE='I';
39 if (N<50){
40 dstev(&JOBZ,&NN,D,E,Z,&LDZ,WORK,&INFO);
41 } else {
42
43 if (N<=1) IL=1;
44 if (N<=1) IU=1;
45
46 LWORK=20*N+1;
47 LIWORK=10*N+1;
48
49 dstevr(&JOBZ,&RANGE,&NN,D,E,&VL,&VU,&IL,&IU,&ABSTOL,&M,WW,Z,&LDZ,ISUPPZ,WORK,&LWORK,IWORK,&LIWORK,&INFO);
50
51 if (N==1){
52 D[0]=WW[0];
53 } else if (N>=2){
54 D[N-2]=WW[0];
55 D[N-1]=WW[1];
56 }
57
58 }
59 return INFO;
60}
61
62/* static int id1=0, id2=0; */
63#undef __FUNCT__
64#define __FUNCT__ "MatMult3"
65static int MatMult3(Mat3 A, SDPConeVec X, SDPConeVec Y){
66
67 int info=0;
68 double minus_one=-1.0;
69
70 DSDPFunctionBegin;
71 /* DSDPEventLogBegin(id2); */
72 if (A->type==2){
73 info=DSDPVMatMult(A->x,X,Y);DSDPCHKERR(info);
74 } else {
75 info=DSDPDualMatCholeskySolveBackward(A->ss,X,Y); DSDPCHKERR(info);
76 info=DSDPDSMatMult(A->ds,Y,A->V); DSDPCHKERR(info);
77 info=DSDPDualMatCholeskySolveForward(A->ss,A->V,Y); DSDPCHKERR(info);
78 info=SDPConeVecScale(minus_one,Y); DSDPCHKERR(info);
79 }
80 /* DSDPEventLogEnd(id2);*/
81 DSDPFunctionReturn(0);
82}
83
84
85#undef __FUNCT__
86#define __FUNCT__ "DSDPLanczosInitialize"
93 DSDPFunctionBegin;
94 /* Build Lanczos structures */
95 LZ->n=0;
96 LZ->lanczosm=0;
97 LZ->maxlanczosm=20;
98 LZ->type=0;
99 LZ->dwork4n=0;
100 LZ->Q=0;
101 LZ->darray=0;
102 /*
103 if (id1==0 && id2==0){
104 DSDPEventLogRegister("STEP EIGS",&id1); printf("ID1: %d\n",id1);
105 DSDPEventLogRegister("STEP MULT",&id2); printf("ID2: %d\n",id2);
106 }
107 */
108 DSDPFunctionReturn(0);
109}
110
111#undef __FUNCT__
112#define __FUNCT__ "DSDPSetMaximumLanczosIterations"
120 DSDPFunctionBegin;
121 if (maxlanczos>0) LZ->lanczosm=maxlanczos;
122 DSDPFunctionReturn(0);
123}
124
125#undef __FUNCT__
126#define __FUNCT__ "DSDPFastLanczosSetup"
134 int i,n,info;
135 DSDPFunctionBegin;
136 /* Build Lanczos structures */
137 info=SDPConeVecGetSize(V,&n);DSDPCHKERR(info);
138 LZ->lanczosm=DSDPMin(LZ->maxlanczosm,n);
139 LZ->type=1;
140 LZ->n=n;
141 if (LZ->lanczosm<50){
142 DSDPCALLOC2(&LZ->dwork4n,double,4*(LZ->lanczosm)+2,&info); DSDPCHKERR(info);
143 DSDPCALLOC2(&LZ->iwork10n,int,1,&info); DSDPCHKERR(info);
144 } else {
145 DSDPCALLOC2(&LZ->dwork4n,double,23*(LZ->lanczosm)+2,&info); DSDPCHKERR(info);
146 DSDPCALLOC2(&LZ->iwork10n,int,10*(LZ->lanczosm),&info); DSDPCHKERR(info);
147 }
148 DSDPCALLOC2(&LZ->Q,SDPConeVec,2,&info);DSDPCHKERR(info);
149 for (i=0;i<2;i++){
150 info = SDPConeVecDuplicate(V,&LZ->Q[i]);DSDPCHKERR(info);
151 }
152 DSDPFunctionReturn(0);
153}
154
155#undef __FUNCT__
156#define __FUNCT__ "DSDPRobustLanczosSetup"
164 int i,n,leig,info;
165 DSDPFunctionBegin;
166 /* Build Lanczos structures */
167 info=SDPConeVecGetSize(V,&n);DSDPCHKERR(info);
168 leig=DSDPMin(LZ->maxlanczosm,n);
169 LZ->n=n;
170 LZ->lanczosm=leig;
171 LZ->type=2;
172
173 DSDPCALLOC2(&LZ->dwork4n,double,3*(LZ->lanczosm)+1,&info); DSDPCHKERR(info);
174 DSDPCALLOC2(&LZ->darray,double,(leig*leig),&info); DSDPCHKERR(info);
175 DSDPCALLOC2(&LZ->Q,SDPConeVec,(leig+1),&info);DSDPCHKERR(info);
176
177 for (i=0;i<=leig;i++){
178 info = SDPConeVecDuplicate(V,&LZ->Q[i]);DSDPCHKERR(info);
179 }
180 info = SDPConeVecCreate(leig,&LZ->Tv);DSDPCHKERR(info);
181 DSDPFunctionReturn(0);
182}
183
184#undef __FUNCT__
185#define __FUNCT__ "DSDPLanczosDestroy"
192 int i,info;
193 DSDPFunctionBegin;
194 if (LZ->type==2){
195 for (i=0;i<=LZ->lanczosm;i++){
196 info = SDPConeVecDestroy(&LZ->Q[i]);DSDPCHKERR(info);
197 }
198 info=SDPConeVecDestroy(&LZ->Tv);DSDPCHKERR(info);
199 DSDPFREE(&LZ->darray,&info);DSDPCHKERR(info);
200 } else if (LZ->type==1){
201 info = SDPConeVecDestroy(&LZ->Q[1]);DSDPCHKERR(info);
202 info = SDPConeVecDestroy(&LZ->Q[0]);DSDPCHKERR(info);
203 DSDPFREE(&LZ->iwork10n,&info);DSDPCHKERR(info);
204 }
205 DSDPFREE(&LZ->Q,&info);DSDPCHKERR(info);
206 DSDPFREE(&LZ->dwork4n,&info);DSDPCHKERR(info);
207 info=DSDPLanczosInitialize(LZ);DSDPCHKERR(info);
208 DSDPFunctionReturn(0);
209}
210
211#undef __FUNCT__
212#define __FUNCT__ "DSDPLanczosMinXEig"
213int DSDPLanczosMinXEig( DSDPLanczosStepLength *LZ, DSDPVMat X, SDPConeVec W1, SDPConeVec W2, double *mineig ){
214 int info,m;
215 double smaxstep;
216 struct _P_Mat3 PP;
217 Mat3 A=&PP;
218
219 DSDPFunctionBegin;
220 A->type=2;
221 A->x=X;
222 A->V=W2;
223 m=LZ->lanczosm;
224
225 if (LZ->type==1){
226 info = ComputeStepFAST(A,LZ->Q,m,W1,LZ->dwork4n,LZ->iwork10n,&smaxstep,mineig);DSDPCHKERR(info);
227 } else if (LZ->type==2){
228 info = ComputeStepROBUST(A,LZ->Q,m,LZ->Q[m],W1,LZ->darray/*LZ->TT*/,LZ->Tv,LZ->dwork4n,&smaxstep,mineig);DSDPCHKERR(info);
229 } else {
230 DSDPSETERR1(1,"Lanczos Step Length Has not been SetUp. Type: %d\n",LZ->type);
231 }
232 DSDPFunctionReturn(0);
233}
234
235#undef __FUNCT__
236#define __FUNCT__ "DSDPLanczosStepSize"
248 int info,m;
249 double smaxstep,mineig;
250 struct _P_Mat3 PP;
251 Mat3 A=&PP;
252
253 DSDPFunctionBegin;
254 A->ss=S;
255 A->ds=DS; A->V=W2;
256 A->type=1;
257 m=LZ->lanczosm;
258
259 if (LZ->type==1){
260 info = ComputeStepFAST(A,LZ->Q,m,W1,LZ->dwork4n,LZ->iwork10n,&smaxstep,&mineig);DSDPCHKERR(info);
261 *maxstep=smaxstep;
262 } else if (LZ->type==2){
263 info = ComputeStepROBUST(A,LZ->Q,m,LZ->Q[m],W1,LZ->darray/*LZ->TT*/,LZ->Tv,LZ->dwork4n,&smaxstep,&mineig);DSDPCHKERR(info);
264 *maxstep=smaxstep;
265 } else {
266 DSDPSETERR1(1,"Lanczos Step Length Has not been SetUp. Type: %d\n",LZ->type);
267 }
268 DSDPFunctionReturn(0);
269}
270
271
272
273#undef __FUNCT__
274#define __FUNCT__ "ComputeStepROBUST"
275static int ComputeStepROBUST(Mat3 A, SDPConeVec *Q, int m, SDPConeVec W, SDPConeVec R, double*darray, SDPConeVec QAQTv, double *dwork, double *maxstep , double *mineig){
276
277 int i,j,n,info;
278 double tt,wnorm, phi;
279 double lambda1=0,lambda2=0,delta=0;
280 double res1,res2,beta;
281 double one=1.0;
282 double *eigvec;
283 int N, LDA, LWORK;
284
285 DSDPFunctionBegin;
286
287 memset((void*)darray,0,m*m*sizeof(double));
288 if (A->type==1){
289 for (i=0; i<m; i++){ darray[i*m+i]=-1.0;}
290 } else {
291 for (i=0; i<m; i++){ darray[i*m+i]=1.0;}
292 }
293
294 info = SDPConeVecSet(one,Q[0]);DSDPCHKERR(info);
295 info = SDPConeVecNormalize(Q[0]);DSDPCHKERR(info);
296
297 for (i=0; i<m; i++){
298 info = MatMult3(A,Q[i],W);DSDPCHKERR(info);
299 info = SDPConeVecNorm2(W,&phi);DSDPCHKERR(info);
300 if (phi!=phi){ *maxstep = 0.0; return 0;}
301 if (i>0){
302 tt=-darray[i*m+i-1];
303 info = SDPConeVecAXPY(tt,Q[i-1],W);DSDPCHKERR(info);
304 }
305 info = SDPConeVecDot(W,Q[i],&tt);DSDPCHKERR(info);
306 darray[i*m+i]=tt;
307 tt*=-1.0;
308 info = SDPConeVecAXPY(tt,Q[i],W);DSDPCHKERR(info);
309 info = SDPConeVecNorm2(W,&wnorm);DSDPCHKERR(info);
310 if (wnorm <= 0.8 * phi){
311 for (j=0;j<=i;j++){
312 info = SDPConeVecDot(W,Q[j],&tt);DSDPCHKERR(info);
313 if (tt==tt){tt*=-1.0;} else {tt=0;}
314 info = SDPConeVecAXPY(tt,Q[j],W);DSDPCHKERR(info);
315 darray[j*m+i]-=tt;
316 if (i!=j){ darray[i*m+j]-=tt; }
317 }
318 }
319
320 info = SDPConeVecNorm2(W,&wnorm);DSDPCHKERR(info);
321 if (i<m-1){
322 darray[i*m+i+1]=wnorm;
323 darray[i*m+m+i]=wnorm;
324 }
325 if (fabs(wnorm)<=1.0e-14) break;
326 info=SDPConeVecCopy(W,Q[i+1]);DSDPCHKERR(info);
327 info=SDPConeVecNormalize(Q[i+1]); DSDPCHKERR(info);
328
329 }
330 /*
331 DSDPLogInfo("Step Length: Lanczos Iterates: %2.0f, ",i);
332 DSDPLogInfo("VNorm: %3.1e, ",wnorm);
333 */
334 /*
335 printf(" --- TRI DIAGONAL MATRIX ---- \n");
336 */
337
338
339 LWORK=DSDPMax(3*m,1); LDA=DSDPMax(1,m); N=m;
340 info=SDPConeVecGetArray(QAQTv,&eigvec);DSDPCHKERR(info);
341 info=DSDPGetEigsSTEP(darray,m,0,0,0,0,eigvec,m,dwork,LWORK,0,0);DSDPCHKERR(info);
342 info=SDPConeVecRestoreArray(QAQTv,&eigvec);DSDPCHKERR(info);
343
344 if (N==0){
345 lambda1=-0.0;
346 delta=1.0e-20;
347 *mineig=0;
348 } else if (N==1){
349 info=SDPConeVecGetArray(QAQTv,&eigvec);DSDPCHKERR(info);
350 lambda1=-eigvec[0];
351 info=SDPConeVecRestoreArray(QAQTv,&eigvec);DSDPCHKERR(info);
352 delta=1.0e-20;
353 *mineig=lambda1;
354 } else if (N>1){
355 info=SDPConeVecGetArray(QAQTv,&eigvec);DSDPCHKERR(info);
356 *mineig=eigvec[0];
357 lambda1=-eigvec[N-1];
358 lambda2=-eigvec[N-2];
359 info=SDPConeVecRestoreArray(QAQTv,&eigvec);DSDPCHKERR(info);
360
361 info = SDPConeVecZero(W);DSDPCHKERR(info);
362 for (i=0;i<m;i++){
363 tt=darray[(N-1)*m+i];
364 info=SDPConeVecAXPY(tt,Q[i],W);DSDPCHKERR(info);
365 }
366 info = MatMult3(A,W,R);DSDPCHKERR(info);
367 info = SDPConeVecAXPY(lambda1,W,R);DSDPCHKERR(info);
368 info = SDPConeVecNorm2(R,&res1);DSDPCHKERR(info);
369
370 info = SDPConeVecZero(W);DSDPCHKERR(info);
371 for (i=0;i<m;i++){
372 tt=darray[(N-2)*m+i];
373 info=SDPConeVecAXPY(tt,Q[i],W);DSDPCHKERR(info);
374 }
375 info = MatMult3(A,W,R);DSDPCHKERR(info);
376 info = SDPConeVecAXPY(lambda2,W,R);DSDPCHKERR(info);
377 info = SDPConeVecNorm2(R,&res2);DSDPCHKERR(info);
378
379 tt = -lambda1 + lambda2 - res2;
380 if (tt>0) beta=tt;
381 else beta=1.0e-20;
382 delta = DSDPMin(res1,sqrt(res1)/beta);
383
384 }
385
386
387 if (delta-lambda1>0)
388 *maxstep = 1.0/(delta-lambda1);
389 else
390 *maxstep = 1.0e+30;
391
392 info=SDPConeVecGetSize(W,&n);DSDPCHKERR(info);
393 DSDPLogInfo(0,19,"Robust Lanczos StepLength: Iterates %d, Max: %d, BlockSize: %d, Lambda1: %4.2e, Res1: %4.2e, Lambda2: %4.2e, Res2: %4.2e, Delta: %4.2e, MaxStep: %4.2e\n",i,m,n,lambda1,res1*res1,lambda2,res2*res2,delta,*maxstep);
394
395 DSDPFunctionReturn(0);
396}
397
398#undef __FUNCT__
399#define __FUNCT__ "ComputeStepFAST"
400static int ComputeStepFAST(Mat3 A, SDPConeVec *Q, int m, SDPConeVec W, double *dwork, int *iwork,double *maxstep ,double *mineig){
401
402 int i,j,n,info;
403 double tt,wnorm, phi;
404 double lambda1=0,lambda2=0,delta=0;
405 double res1,res2,beta;
406 double one=1.0;
407 int N=m;
408 double *diag,*subdiag,*ddwork;
409
410 DSDPFunctionBegin;
411 diag=dwork;
412 subdiag=dwork+m;
413 ddwork=dwork+2*m;
414
415 if (A->type==1){
416 for (i=0; i<m; i++){ diag[i]=-1; subdiag[i]=0;}
417 } else {
418 for (i=0; i<m; i++){ diag[i]=1.0; subdiag[i]=0;}
419 }
420 info = SDPConeVecSet(one,Q[0]);DSDPCHKERR(info);
421 info = SDPConeVecNormalize(Q[0]);DSDPCHKERR(info);
422
423 for (i=0; i<m; i++){
424 info = MatMult3(A,Q[0],W);DSDPCHKERR(info);
425 info = SDPConeVecNorm2(W,&phi);DSDPCHKERR(info);
426 if (phi!=phi){ *maxstep = 0.0; return 0;}
427 if (i>0){
428 tt=-subdiag[i-1];
429 info = SDPConeVecAXPY(tt,Q[1],W);DSDPCHKERR(info);
430 }
431 info = SDPConeVecDot(W,Q[0],&tt);DSDPCHKERR(info);
432 diag[i]=tt;
433 tt*=-1.0;
434 info = SDPConeVecAXPY(tt,Q[0],W);DSDPCHKERR(info);
435 info = SDPConeVecNorm2(W,&wnorm);DSDPCHKERR(info);
436 if (wnorm <= 1.0 * phi){
437 for (j=0;j<=i;j++){
438 if (j==i-1){
439 info = SDPConeVecDot(W,Q[1],&tt);DSDPCHKERR(info);
440 if (tt==tt){tt*=-1.0;} else {tt=0;}
441 info = SDPConeVecAXPY(tt,Q[1],W);DSDPCHKERR(info);
442 subdiag[i-1]-=tt;
443 } else if (j==i){
444 info = SDPConeVecDot(W,Q[0],&tt);DSDPCHKERR(info);
445 if (tt==tt){tt*=-1.0;} else {tt=0;}
446 info = SDPConeVecAXPY(tt,Q[0],W);DSDPCHKERR(info);
447 diag[i]-=tt;
448 }
449
450 }
451 }
452
453 info = SDPConeVecNorm2(W,&wnorm);DSDPCHKERR(info);
454 /* printf("PHI: %4.4e, VNORM: %4.2e Diag: %4.2e\n",phi,wnorm,diag[i]); */
455 if (i<m-1){
456 subdiag[i]=wnorm;
457 }
458 if (fabs(wnorm)<=1.0e-10){i++;break;}
459 info=SDPConeVecCopy(Q[0],Q[1]);DSDPCHKERR(info);
460 info=SDPConeVecCopy(W,Q[0]);DSDPCHKERR(info);
461 info=SDPConeVecNormalize(Q[0]); DSDPCHKERR(info);
462
463 }
464
465 /* DSDPEventLogBegin(id1); */
466 info=DSDPGetTriDiagonalEigs(m,diag,subdiag,ddwork,iwork); DSDPCHKERR(info);
467 /* DSDPEventLogEnd(id1); */
468 if (N==0){
469 lambda1=-0.0;
470 delta=1.0e-20;
471 *mineig=0;
472 } else if (N==1){
473 lambda1=-diag[0];
474 delta=1.0e-20;
475 *mineig=diag[0];
476 } else if (N>1){
477 lambda1=-diag[N-1];
478 lambda2=-diag[N-2];
479
480 res1=1.0e-8;
481 res2=1.0e-8;
482
483 tt = -lambda1 + lambda2 - res2;
484 if (tt>0) beta=tt;
485 else beta=1.0e-20;
486 delta = DSDPMin(res1,sqrt(res1)/beta);
487
488 *mineig=diag[0];
489 }
490
491
492 if (delta-lambda1>0)
493 *maxstep = 1.0/(delta-lambda1);
494 else
495 *maxstep = 1.0e+30;
496
497 info=SDPConeVecGetSize(W,&n);DSDPCHKERR(info);
498 DSDPLogInfo(0,19,"Step Length: Fast Lanczos Iterates: %2d, Max: %d, Block Size: %d, VNorm: %3.1e, Lambda1: %4.4e, Lambda2: %4.4e, Delta: %4.2e, Maxstep: %4.2e\n",
499 i,m,n,wnorm,lambda1,lambda2,delta,*maxstep);
500
501 DSDPFunctionReturn(0);
502}
int DSDPDSMatMult(DSDPDSMat A, SDPConeVec X, SDPConeVec Y)
Set values into the matrix.
Definition: dsdpdsmat.c:154
The interface between the SDPCone and the Delta S matrix.
int DSDPDualMatCholeskySolveBackward(DSDPDualMat S, SDPConeVec B, SDPConeVec X)
Backward triangular solve.
Definition: dsdpdualmat.c:295
int DSDPDualMatCholeskySolveForward(DSDPDualMat S, SDPConeVec B, SDPConeVec X)
Forward triangular solve.
Definition: dsdpdualmat.c:267
The interface between the SDPCone and the matrix S.
Lanczos procedure determines the maximum step length.
DSDP uses BLAS and LAPACK for many of its operations.
int DSDPSetMaximumLanczosIterations(DSDPLanczosStepLength *LZ, int maxlanczos)
Set parameter.
Definition: dsdpstep.c:119
int DSDPLanczosInitialize(DSDPLanczosStepLength *LZ)
Initialize Lanczos structure.
Definition: dsdpstep.c:92
int DSDPLanczosDestroy(DSDPLanczosStepLength *LZ)
Free data structure.
Definition: dsdpstep.c:191
int DSDPRobustLanczosSetup(DSDPLanczosStepLength *LZ, SDPConeVec V)
Use slowerer but more robust method.
Definition: dsdpstep.c:163
int DSDPFastLanczosSetup(DSDPLanczosStepLength *LZ, SDPConeVec V)
Use Lanczos procedure. Assume off tridiagonal entries are zero.
Definition: dsdpstep.c:133
int DSDPLanczosStepSize(DSDPLanczosStepLength *LZ, SDPConeVec W1, SDPConeVec W2, DSDPDualMat S, DSDPDSMat DS, double *maxstep)
Compute distance to boundary.
Definition: dsdpstep.c:247
Error handling, printing, and profiling.
int DSDPVMatMult(DSDPVMat X, SDPConeVec Z, SDPConeVec Y)
Multiply X by a vector.
Definition: dsdpxmat.c:301
The interface between the SDPCone and the dense matrix array.
int SDPConeVecAXPY(double alpha, SDPConeVec x, SDPConeVec y)
Add a multiple of X to Y.
Definition: sdpconevec.c:178
int SDPConeVecDuplicate(SDPConeVec V1, SDPConeVec *V2)
Allocate another vector with the same structure as the first.
Definition: sdpconevec.c:195
int SDPConeVecNorm2(SDPConeVec VV, double *vnorm)
Compute the Euclidean norm.
Definition: sdpconevec.c:143
int SDPConeVecCopy(SDPConeVec v1, SDPConeVec v2)
Copy v1 to v2.
Definition: sdpconevec.c:103
int SDPConeVecScale(double alpha, SDPConeVec VV)
Compute the Euclidean norm.
Definition: sdpconevec.c:161
int SDPConeVecZero(SDPConeVec V)
Zero the elements of the vector.
Definition: sdpconevec.c:67
int SDPConeVecDot(SDPConeVec V1, SDPConeVec V2, double *ans)
Inner product of two vectors.
Definition: sdpconevec.c:125
int SDPConeVecNormalize(SDPConeVec V)
Scale the vector to norm of 1.
Definition: sdpconevec.c:84
int SDPConeVecSet(double alpha, SDPConeVec V)
Set each element of vector to this number.
Definition: sdpconevec.c:211
Symmetric Delta S matrix for one block in the semidefinite cone.
Definition: dsdpdsmat.h:23
Represents an S matrix for one block in the semidefinite cone.
Definition: dsdpdualmat.h:18
Apply Lanczos prodedure to find distance to boundary.
Definition: dsdplanczos.h:13
Dense symmetric matrix for one block in the semidefinite cone.
Definition: dsdpxmat.h:17
Vector whose length corresponds to dimension of a block in a cone.
Definition: sdpconevec.h:13