DSDP
dsdpcops.c
Go to the documentation of this file.
1#include "dsdp.h"
2#include "dsdpsys.h"
8#define DSDPCHKCONEERR(a,b); { if (b){ DSDPSETERR1(b,"Cone Number: %d,\n",a);} }
9
10static int ConeSetup=0,ConeComputeS=0,ConeComputeSS=0,ConeComputeH=0,ConeHMultiplyAdd=0,ConeMaxPStep=0,ConeMaxDStep=0,ConePotential=0,ConeComputeX=0,ConeView=0,ConeDestroy=0,ConeXEigs=0,ConeRHS=0,ConeInvertS=0;
11static int DSDPRegisterConeEvents(void);
12int DSDPSetSchurMatOps(DSDP,struct DSDPSchurMat_Ops*, void*);
14/*
15int DSDPIdentifySchurColumns(DSDP,int, int*, int*, int);
16*/
17
18#undef __FUNCT__
19#define __FUNCT__ "DSDPZeroConeEvents"
20static int DSDPZeroConeEvents(){
21 DSDPFunctionBegin;
22 ConeSetup=0;ConeComputeS=0;ConeComputeSS=0;ConeComputeH=0;ConeHMultiplyAdd=0;ConeMaxPStep=0;ConeMaxDStep=0;ConePotential=0;ConeComputeX=0;ConeView=0;ConeDestroy=0;ConeXEigs=0;ConeRHS=0;ConeInvertS=0;
23 DSDPFunctionReturn(0);
24}
25
26#undef __FUNCT__
27#define __FUNCT__ "DSDPRegisterConeEvents"
28static int DSDPRegisterConeEvents(){
29
30 DSDPFunctionBegin;
31 DSDPEventLogRegister("Cone Setup 1&2",&ConeSetup);
32 DSDPEventLogRegister("Cone Invert S",&ConeInvertS);
33 DSDPEventLogRegister("Cone RHS",&ConeRHS);
34 DSDPEventLogRegister("Cone Compute Newton Eq.",&ConeComputeH);
35 DSDPEventLogRegister("Cone Newton Multiply-Add",&ConeHMultiplyAdd);
36 DSDPEventLogRegister("Cone Max P Step Length",&ConeMaxPStep);
37 DSDPEventLogRegister("Cone Compute and Factor SP",&ConeComputeSS);
38 DSDPEventLogRegister("Cone Max D Step Length",&ConeMaxDStep);
39 DSDPEventLogRegister("Cone Compute and Factor S",&ConeComputeS);
40 DSDPEventLogRegister("Cone Potential",&ConePotential);
41 DSDPEventLogRegister("Cone View",&ConeView);
42 DSDPEventLogRegister("Cone Compute X",&ConeComputeX);
43 DSDPEventLogRegister("Cone X Residuals",&ConeXEigs);
44 DSDPEventLogRegister("Cone Destroy",&ConeDestroy);
45
46 DSDPFunctionReturn(0);
47}
48
56#undef __FUNCT__
57#define __FUNCT__ "DSDPSetUpCones"
59 int info,kk;
60 DSDPVec yy0=dsdp->y;
61 DSDPFunctionBegin;
62 info=DSDPRegisterConeEvents();
63 DSDPEventLogBegin(ConeSetup);
64 for (kk=0;kk<dsdp->ncones;kk++){
65 DSDPEventLogBegin(dsdp->K[kk].coneid);
66 info=DSDPConeSetUp(dsdp->K[kk].cone,yy0);DSDPCHKCONEERR(kk,info);
67 DSDPEventLogEnd(dsdp->K[kk].coneid);
68 }
69 DSDPEventLogEnd(ConeSetup);
70 DSDPFunctionReturn(0);
71}
72
82#undef __FUNCT__
83#define __FUNCT__ "DSDPSetUpCones2"
85 int info,kk;
86 DSDPFunctionBegin;
87 DSDPEventLogBegin(ConeSetup);
88 for (kk=0;kk<dsdp->ncones;kk++){
89 DSDPEventLogBegin(dsdp->K[kk].coneid);
90 info=DSDPConeSetUp2(dsdp->K[kk].cone,yy0,M);DSDPCHKCONEERR(kk,info);
91 DSDPEventLogEnd(dsdp->K[kk].coneid);
92 }
93 DSDPEventLogEnd(ConeSetup);
94 DSDPFunctionReturn(0);
95}
96
97
105#undef __FUNCT__
106#define __FUNCT__ "DSDPDestroyCones"
108 int info,kk,ncones=dsdp->ncones;
109 DSDPFunctionBegin;
110 DSDPEventLogBegin(ConeDestroy);
111 for (kk=ncones-1;kk>=0; kk--){
112 DSDPEventLogBegin(dsdp->K[kk].coneid);
113 info=DSDPConeDestroy(&dsdp->K[kk].cone);DSDPCHKCONEERR(kk,info);
114 DSDPEventLogEnd(dsdp->K[kk].coneid);
115 info=DSDPConeInitialize(&dsdp->K[kk].cone);DSDPCHKCONEERR(kk,info);
116 dsdp->ncones--;
117 }
118 if (dsdp->maxcones>0){
119 DSDPFREE(&dsdp->K,&info);DSDPCHKERR(info);
120 dsdp->K=0;
121 dsdp->maxcones=0;
122 }
123 DSDPEventLogEnd(ConeDestroy);
124 info=DSDPZeroConeEvents();DSDPCHKERR(info);
125 DSDPFunctionReturn(0);
126}
127
128
140#undef __FUNCT__
141#define __FUNCT__ "DSDPComputeHessian"
143 int info,kk; double r;
144 DSDPFunctionBegin;
145 DSDPEventLogBegin(ConeComputeH);
146 dsdp->schurmu=dsdp->mutarget;
147 info=DSDPVecGetR(dsdp->y,&r);DSDPCHKERR(info);
148 info=DSDPSchurMatSetR(dsdp->M,r);DSDPCHKERR(info);
149 info=DSDPSchurMatZeroEntries(M);DSDPCHKERR(info);
150 info=DSDPVecZero(vrhs1);DSDPCHKERR(info);
151 info=DSDPVecZero(vrhs2);DSDPCHKERR(info);
152 info=DSDPVecZero(M.schur->rhs3);DSDPCHKERR(info);
153 info=DSDPObjectiveGH(dsdp,M,vrhs1); DSDPCHKERR(info);
154 for (kk=dsdp->ncones-1;kk>=0;kk--){
155 DSDPEventLogBegin(dsdp->K[kk].coneid);
156 info=DSDPConeComputeHessian(dsdp->K[kk].cone,dsdp->schurmu,M,vrhs1,vrhs2);DSDPCHKCONEERR(kk,info);
157 DSDPEventLogEnd(dsdp->K[kk].coneid);
158 }
159 info=DSDPSchurMatAssemble(M);DSDPCHKERR(info);
160 /* DSDPSchurMatView(M); */
161 info=DSDPSchurMatReducePVec(M,vrhs1);DSDPCHKERR(info);
162 info=DSDPSchurMatReducePVec(M,vrhs2);DSDPCHKERR(info);
163 info=DSDPSchurMatReducePVec(M,M.schur->rhs3);DSDPCHKERR(info);
164 if (0 && dsdp->UsePenalty==DSDPNever){
165 info=DSDPVecAXPY(1.0,M.schur->rhs3,vrhs2);DSDPCHKERR(info);
166 info=DSDPVecZero(M.schur->rhs3);DSDPCHKERR(info);
167 info=DSDPVecZero(M.schur->dy3);DSDPCHKERR(info);
168 info=DSDPVecSetR(vrhs1,0);DSDPCHKERR(info);
169 info=DSDPVecSetR(vrhs2,r);DSDPCHKERR(info);
170 }
171 DSDPEventLogEnd(ConeComputeH);
172 DSDPFunctionReturn(0);
173}
174
175
176#undef __FUNCT__
177#define __FUNCT__ "DSDPHessianMultiplyAdd"
189 int info,kk;
190 DSDPVec vrow=dsdp->sles->BR;
191 DSDPFunctionBegin;
192 DSDPEventLogBegin(ConeHMultiplyAdd);
193
194 info=DSDPSchurMatRowScaling(dsdp->M,vrow);DSDPCHKERR(info);
195 for (kk=0;kk<dsdp->ncones;kk++){
196 DSDPEventLogBegin(dsdp->K[kk].coneid);
197 info=DSDPConeMultiplyAdd(dsdp->K[kk].cone,dsdp->schurmu,vrow,v,vv);DSDPCHKCONEERR(kk,info);
198 DSDPEventLogEnd(dsdp->K[kk].coneid);
199 }
200 info=DSDPSchurMatReducePVec(dsdp->M,vv);DSDPCHKERR(info);
201 DSDPEventLogEnd(ConeHMultiplyAdd);
202 DSDPFunctionReturn(0);
203}
204
205#undef __FUNCT__
206#define __FUNCT__ "DSDPComputeG"
215int DSDPComputeG( DSDP dsdp , DSDPVec vt, DSDPVec vrhs1, DSDPVec vrhs2){
216 int info,kk; double r;
217 DSDPFunctionBegin;
218 DSDPEventLogBegin(ConeRHS);
219 info=DSDPVecZero(vrhs1);DSDPCHKERR(info);
220 info=DSDPVecZero(vrhs2);DSDPCHKERR(info);
221 info=DSDPVecGetR(dsdp->y,&r);DSDPCHKERR(info);
222 info=DSDPSchurMatSetR(dsdp->M,r);DSDPCHKERR(info);
223 info=DSDPSchurMatRowScaling(dsdp->M,vt);DSDPCHKERR(info);
224 info=DSDPObjectiveGH(dsdp,dsdp->M,vrhs1); DSDPCHKERR(info);
225 if (0 && r==0){info=DSDPVecSetR(vrhs1,0);info=DSDPVecSetR(vrhs2,0);}
226 /* info=DSDPVecScale(1.0/dsdp->schurmu,vrhs1); DSDPCHKERR(info); */
227 for (kk=0;kk<dsdp->ncones;kk++){
228 DSDPEventLogBegin(dsdp->K[kk].coneid);
229 info=DSDPConeComputeRHS(dsdp->K[kk].cone,dsdp->schurmu,vt,vrhs1,vrhs2);DSDPCHKCONEERR(kk,info);
230 DSDPEventLogEnd(dsdp->K[kk].coneid);
231 }
232 DSDPEventLogEnd(ConeRHS);
233 info=DSDPSchurMatReducePVec(dsdp->M,vrhs1);DSDPCHKERR(info);
234 info=DSDPSchurMatReducePVec(dsdp->M,vrhs2);DSDPCHKERR(info);
235 DSDPFunctionReturn(0);
236}
237
238#undef __FUNCT__
239#define __FUNCT__ "DSDPComputeANorm2"
246int DSDPComputeANorm2( DSDP dsdp , DSDPVec Anorm2){
247 int info,kk;
248 DSDPFunctionBegin;
249 for (kk=0;kk<dsdp->ncones;kk++){
250 DSDPEventLogBegin(dsdp->K[kk].coneid);
251 info=DSDPConeANorm2(dsdp->K[kk].cone,Anorm2);DSDPCHKCONEERR(kk,info);
252 DSDPEventLogEnd(dsdp->K[kk].coneid);
253 }
254 DSDPFunctionReturn(0);
255}
256
257
270#undef __FUNCT__
271#define __FUNCT__ "DSDPComputeSS"
272int DSDPComputeSS(DSDP dsdp, DSDPVec Y, DSDPDualFactorMatrix flag, DSDPTruth *ispsdefinite){
273 int info,kk;
275 DSDPFunctionBegin;
276 if (flag==DUAL_FACTOR){
277 DSDPEventLogBegin(ConeComputeS);
278 } else if (flag==PRIMAL_FACTOR){
279 DSDPEventLogBegin(ConeComputeSS);
280 }
281 for (kk=dsdp->ncones-1; kk>=0 && psd==DSDP_TRUE;kk--){
282 DSDPEventLogBegin(dsdp->K[kk].coneid);
283 info=DSDPConeComputeS(dsdp->K[kk].cone,Y,flag,&psd); DSDPCHKCONEERR(kk,info);
284 DSDPEventLogEnd(dsdp->K[kk].coneid);
285 }
286 *ispsdefinite=psd;
287 if (flag==DUAL_FACTOR){
288 DSDPEventLogEnd(ConeComputeS);
289 } else if (flag==PRIMAL_FACTOR){
290 DSDPEventLogEnd(ConeComputeSS);
291 }
292 DSDPFunctionReturn(0);
293}
294
305#undef __FUNCT__
306#define __FUNCT__ "DSDPInvertS"
308 int info,kk;
309 DSDPFunctionBegin;
310 DSDPEventLogBegin(ConeInvertS);
311 for (kk=0;kk<dsdp->ncones;kk++){
312 DSDPEventLogBegin(dsdp->K[kk].coneid);
313 info=DSDPConeInvertS(dsdp->K[kk].cone); DSDPCHKCONEERR(kk,info);
314 DSDPEventLogEnd(dsdp->K[kk].coneid);
315 }
316 DSDPEventLogEnd(ConeInvertS);
317 DSDPFunctionReturn(0);
318}
319
334#undef __FUNCT__
335#define __FUNCT__ "DSDPComputeMaxStepLength"
336int DSDPComputeMaxStepLength(DSDP dsdp, DSDPVec DY, DSDPDualFactorMatrix flag, double *maxsteplength){
337 int info,kk;
338 double msteplength=1.0e30,conesteplength;
339 DSDPFunctionBegin;
340 if (flag==DUAL_FACTOR){
341 DSDPEventLogBegin(ConeMaxDStep);
342 } else if (flag==PRIMAL_FACTOR){
343 DSDPEventLogBegin(ConeMaxPStep);
344 }
345 for (kk=0;kk<dsdp->ncones;kk++){
346 DSDPEventLogBegin(dsdp->K[kk].coneid);
347 conesteplength=1.0e20;
348 info=DSDPConeComputeMaxStepLength(dsdp->K[kk].cone,DY,flag,&conesteplength);DSDPCHKCONEERR(kk,info);
349 msteplength=DSDPMin(msteplength,conesteplength);
350 DSDPEventLogEnd(dsdp->K[kk].coneid);
351 }
352 *maxsteplength=msteplength;
353 if (flag==DUAL_FACTOR){
354 DSDPEventLogEnd(ConeMaxDStep);
355 } else if (flag==PRIMAL_FACTOR){
356 DSDPEventLogEnd(ConeMaxPStep);
357 }
358 DSDPFunctionReturn(0);
359}
360
361
376#undef __FUNCT__
377#define __FUNCT__ "DSDPPassXVectors"
378int DSDPPassXVectors(DSDP dsdp, double mu, DSDPVec Y, DSDPVec DY){
379 int info,kk;
380 DSDPFunctionBegin;
381 for (kk=0;kk<dsdp->ncones;kk++){
382 DSDPEventLogBegin(dsdp->K[kk].coneid);
383 info=DSDPConeSetXMaker(dsdp->K[kk].cone,mu,Y,DY);DSDPCHKCONEERR(kk,info);
384 DSDPEventLogEnd(dsdp->K[kk].coneid);
385 }
386 DSDPFunctionReturn(0);
387}
388
389
399#undef __FUNCT__
400#define __FUNCT__ "DSDPGetConicDimension"
401int DSDPGetConicDimension(DSDP dsdp, double *n){
402 int info,kk;
403 double nn,nnn=0;
404 DSDPFunctionBegin;
405 for (kk=0;kk<dsdp->ncones;kk++){
406 nn=0;
407 info=DSDPConeGetDimension(dsdp->K[kk].cone,&nn);DSDPCHKCONEERR(kk,info);
408 nnn+=nn;
409 }
410 *n=nnn;
411 DSDPFunctionReturn(0);
412}
413
414
422#undef __FUNCT__
423#define __FUNCT__ "DSDPViewCones"
425 int info,kk;
426 DSDPFunctionBegin;
427 DSDPEventLogBegin(ConeView);
428 for (kk=0;kk<dsdp->ncones;kk++){
429 DSDPEventLogBegin(dsdp->K[kk].coneid);
430 info=DSDPConeView(dsdp->K[kk].cone);DSDPCHKCONEERR(kk,info);
431 DSDPEventLogEnd(dsdp->K[kk].coneid);
432 }
433 DSDPEventLogEnd(ConeView);
434 DSDPFunctionReturn(0);
435}
436
437
448#undef __FUNCT__
449#define __FUNCT__ "DSDPMonitorCones"
450int DSDPMonitorCones(DSDP dsdp,int tag){
451 int info,kk;
452 DSDPFunctionBegin;
453 DSDPEventLogBegin(ConeView);
454 for (kk=0;kk<dsdp->ncones;kk++){
455 DSDPEventLogBegin(dsdp->K[kk].coneid);
456 info=DSDPConeMonitor(dsdp->K[kk].cone,tag);DSDPCHKCONEERR(kk,info);
457 DSDPEventLogEnd(dsdp->K[kk].coneid);
458 }
459 DSDPEventLogEnd(ConeView);
460 DSDPFunctionReturn(0);
461}
462
463
472#undef __FUNCT__
473#define __FUNCT__ "DSDPSparsityInSchurMat"
474int DSDPSchurSparsity(DSDP dsdp, int row, int rnnz[], int m){
475 int info,kk;
476 DSDPFunctionBegin;
477 for (kk=0;kk<dsdp->ncones;kk++){
478 /* DSDPEventLogBegin(dsdp->K[kk].coneid); */
479 info=DSDPConeSparsityInSchurMat(dsdp->K[kk].cone,row,rnnz,m+2);DSDPCHKCONEERR(kk,info);
480 /* DSDPEventLogEnd(dsdp->K[kk].coneid); */
481 }
482 DSDPFunctionReturn(0);
483}
484
493#undef __FUNCT__
494#define __FUNCT__ "DSDPComputeLogSDeterminant"
495int DSDPComputeLogSDeterminant(DSDP dsdp, double *logdet){
496 int info,kk;
497 double coneobjective,conepotential,llogdet=0;
498 DSDPFunctionBegin;
499 DSDPEventLogBegin(ConePotential);
500 for (kk=0;kk<dsdp->ncones;kk++){
501 DSDPEventLogBegin(dsdp->K[kk].coneid);
502 coneobjective=0;conepotential=0;
503 info=DSDPConeComputeLogSDeterminant(dsdp->K[kk].cone,&coneobjective,&conepotential);DSDPCHKCONEERR(kk,info);
504 llogdet+=conepotential;
505 DSDPEventLogEnd(dsdp->K[kk].coneid);
506 }
507 *logdet=llogdet;
508 DSDPEventLogEnd(ConePotential);
509 DSDPFunctionReturn(0);
510}
511
512
520#undef __FUNCT__
521#define __FUNCT__ "DSDPSetCone"
522int DSDPSetCone(DSDP dsdp,DSDPCone tcone){
523 int info,i,tc;
524 char conename[100];
525 DCone *ccones;
526 DSDPFunctionBegin;
527 if (dsdp->ncones>=dsdp->maxcones){
528 tc=2*dsdp->maxcones+4;
529
530 DSDPCALLOC2(&ccones,DCone,tc,&info);DSDPCHKERR(info);
531 for (i=0;i<dsdp->ncones;i++){ccones[i].cone=dsdp->K[i].cone;}
532 for (i=0;i<dsdp->ncones;i++){ccones[i].coneid=dsdp->K[i].coneid;}
533 DSDPFREE(&dsdp->K,&info);DSDPCHKERR(info);
534 dsdp->K=ccones;
535 dsdp->maxcones=tc;
536 }
537 info=DSDPGetConeName(tcone,conename,100);DSDPCHKERR(info);
538 DSDPEventLogRegister(conename,&tc);
539 dsdp->K[dsdp->ncones].cone=tcone;
540 dsdp->K[dsdp->ncones].coneid=tc;
541 dsdp->ncones++;
542 DSDPFunctionReturn(0);
543}
544
567#undef __FUNCT__
568#define __FUNCT__ "DSDPAddCone"
569int DSDPAddCone(DSDP dsdp,struct DSDPCone_Ops* dsdpops, void* dsdpcone){
570 int info;
571 DSDPCone K;
572 DSDPFunctionBegin;
573 info=DSDPConeInitialize(&K); DSDPCHKERR(info);
574 info=DSDPConeSetData(&K,dsdpops,dsdpcone); DSDPCHKERR(info);
575 info=DSDPSetCone(dsdp,K); DSDPCHKERR(info);
576 DSDPFunctionReturn(0);
577}
578
579
600#undef __FUNCT__
601#define __FUNCT__ "DSDPSetSchurMatOps"
602int DSDPSetSchurMatOps(DSDP dsdp,struct DSDPSchurMat_Ops* sops, void* mdata){
603 int info;
604 DSDPFunctionBegin;
605 info=DSDPSchurMatSetData(&dsdp->M,sops,mdata);DSDPCHKERR(info);
606 DSDPFunctionReturn(0);
607}
608
609
620#undef __FUNCT__
621#define __FUNCT__ "DSDPSetSchurRow"
622int DSDPAddSchurRow(DSDP dsdp,int row, DSDPVec R){
623 int info;
624 DSDPFunctionBegin;
625 info=DSDPSchurMatAddRow(dsdp->M,row,1.0,R);DSDPCHKERR(info);
626 DSDPFunctionReturn(0);
627}
628/*
629#undef __FUNCT__
630#define __FUNCT__ "DSDPIdentifySchurColumns"
631int DSDPIdentifySchurColumns(DSDP dsdp,int row, int *mcol, int *ncols, int m){
632 DSDPFunctionBegin;
633 int info;
634 DSDPVec V;
635 info=DSDPSchurMatRowColumnScaling(dsdp->M,row,V,ncols); DSDPCHKERR(info);
636 DSDPFunctionReturn(1);
637}
638*/
639
652#undef __FUNCT__
653#define __FUNCT__ "DSDPComputeXVariables"
654int DSDPComputeXVariables(DSDP dsdp, double xmakermu, DSDPVec xmakery, DSDPVec xmakerdy, DSDPVec AX, double *tracexs){
655 int kk,info;
656 double ttracexs=0,tttracexs=0,tracex;
657
658 DSDPFunctionBegin;
659 DSDPEventLogBegin(ConeComputeX);
660 info=DSDPVecZero(AX);DSDPCHKERR(info);
661 for (kk=0;kk<dsdp->ncones;kk++){
662 DSDPEventLogBegin(dsdp->K[kk].coneid);
663 tttracexs=0;
664 info=DSDPConeComputeX(dsdp->K[kk].cone,xmakermu,xmakery,xmakerdy,AX,&tttracexs);DSDPCHKCONEERR(kk,info);
665 ttracexs+=tttracexs;
666 DSDPEventLogEnd(dsdp->K[kk].coneid);
667 }
668 info=DSDPVecGetR(AX,&tracex); DSDPCHKERR(info);
669 DSDPLogInfo(0,2,"Trace(X): %4.2e\n",dsdp->tracex);
670 info=DSDPVecAXPY(-1.0,dsdp->b,AX); DSDPCHKERR(info);
671 info=DSDPComputeFixedYX(dsdp->M,AX); DSDPCHKERR(info);
672 *tracexs=ttracexs;
673 info=DSDPVecSetR(AX,tracex); DSDPCHKERR(info);
674 DSDPEventLogEnd(ConeComputeX);
675 DSDPFunctionReturn(0);
676}
677
678
Internal data structure for the DSDP solver.
int DSDPObjectiveGH(DSDP, DSDPSchurMat, DSDPVec)
Compute gradient of dual objective.
Definition: dualimpl.c:381
DSDPDualFactorMatrix
DSDP requires two instances of the data structures S.
@ PRIMAL_FACTOR
@ DUAL_FACTOR
DSDPTruth
Boolean variables.
@ DSDP_TRUE
int DSDPConeView(DSDPCone K)
View contents of the cone.
Definition: dsdpcone.c:358
int DSDPConeComputeX(DSDPCone K, double mu, DSDPVec y, DSDPVec dy, DSDPVec AX, double *tracexs)
Given y,dy, and mu, construct X and add its inner product with the data and S.
Definition: dsdpcone.c:216
int DSDPConeInitialize(DSDPCone *K)
Initialize the pointers to 0.
Definition: dsdpcone.c:495
int DSDPGetConeName(DSDPCone K, char *cname, int maxlength)
Get name of the cone.
Definition: dsdpcone.c:427
int DSDPConeComputeMaxStepLength(DSDPCone K, DSDPVec DY, DSDPDualFactorMatrix flag, double *maxsteplength)
Determine distance to the edge of the cone.
Definition: dsdpcone.c:288
int DSDPConeComputeS(DSDPCone K, DSDPVec Y, DSDPDualFactorMatrix flag, DSDPTruth *ispsdefinite)
Given y, compute S and determine whether its in the cone.
Definition: dsdpcone.c:242
int DSDPConeSetXMaker(DSDPCone K, double mu, DSDPVec y, DSDPVec dy)
Pass information needed to construct X.
Definition: dsdpcone.c:191
int DSDPConeMonitor(DSDPCone K, int tag)
Do anything at in the cone at each iteration.
Definition: dsdpcone.c:380
int DSDPConeComputeLogSDeterminant(DSDPCone K, double *logdetobj, double *logdet)
Evaluate logrithmic barrier function.
Definition: dsdpcone.c:403
int DSDPConeInvertS(DSDPCone K)
Invert the dual matrix S.
Definition: dsdpcone.c:265
int DSDPConeGetDimension(DSDPCone K, double *n)
Provide the dimension of the cone.
Definition: dsdpcone.c:312
int DSDPConeSetUp(DSDPCone K, DSDPVec y)
Factor the data and allocate data structures.
Definition: dsdpcone.c:22
int DSDPConeDestroy(DSDPCone *K)
Free the internal memory of the cone.
Definition: dsdpcone.c:64
int DSDPConeComputeRHS(DSDPCone K, double mu, DSDPVec vrow, DSDPVec rhs1, DSDPVec rhs2)
Compute gradient of barrier function.
Definition: dsdpcone.c:147
int DSDPConeSparsityInSchurMat(DSDPCone K, int row, int rnnz[], int m)
Identify sparsity pattern in a row of the Hessian term.
Definition: dsdpcone.c:338
int DSDPConeComputeHessian(DSDPCone K, double mu, DSDPSchurMat M, DSDPVec vrhs1, DSDPVec vrhs2)
Compute Hessian and gradient of barrier function.
Definition: dsdpcone.c:92
int DSDPConeSetData(DSDPCone *K, struct DSDPCone_Ops *ops, void *data)
Initialize the pointers to 0.
Definition: dsdpcone.c:477
int DSDPConeMultiplyAdd(DSDPCone K, double mu, DSDPVec vrow, DSDPVec v, DSDPVec vv)
Multiply Hessian by a vector and add the result.
Definition: dsdpcone.c:119
int DSDPConeSetUp2(DSDPCone K, DSDPVec yy0, DSDPSchurMat M)
Factor the data and allocate data structures.
Definition: dsdpcone.c:43
int DSDPConeANorm2(DSDPCone K, DSDPVec anorm2)
Add square of 2-norm of data correponding to each variable y.
Definition: dsdpcone.c:168
struct DSDPCone_C DSDPCone
This object holds the data of a SDP, LP, or other cone. Its structure is opaque to the DSDP Solver,...
Definition: dsdpcone.h:27
int DSDPDestroyCones(DSDP dsdp)
Each cone shoudl free its data structures.
Definition: dsdpcops.c:107
int DSDPSetCone(DSDP dsdp, DSDPCone tcone)
Pass a cone to the DSDP solver.
Definition: dsdpcops.c:522
int DSDPViewCones(DSDP dsdp)
Each cone should print its state.
Definition: dsdpcops.c:424
int DSDPComputeLogSDeterminant(DSDP dsdp, double *logdet)
Compute the logarithmic barrier function for the dual varialbe S.
Definition: dsdpcops.c:495
int DSDPSchurSparsity(DSDP dsdp, int row, int rnnz[], int m)
Each cone should print its state.
Definition: dsdpcops.c:474
int DSDPAddCone(DSDP dsdp, struct DSDPCone_Ops *dsdpops, void *dsdpcone)
Apply DSDP to a conic structure.
Definition: dsdpcops.c:569
int DSDPComputeMaxStepLength(DSDP dsdp, DSDPVec DY, DSDPDualFactorMatrix flag, double *maxsteplength)
Compute the maximum step length for the given step direction.
Definition: dsdpcops.c:336
int DSDPSetUpCones(DSDP dsdp)
Each cone should factor data or allocate internal data structures.
Definition: dsdpcops.c:58
int DSDPPassXVectors(DSDP dsdp, double mu, DSDPVec Y, DSDPVec DY)
Pass the information needed to compute the variables X in each cone but do not compute X.
Definition: dsdpcops.c:378
int DSDPComputeHessian(DSDP dsdp, DSDPSchurMat M, DSDPVec vrhs1, DSDPVec vrhs2)
Compute the Schur complement, or Gram, matrix for each cone.
Definition: dsdpcops.c:142
int DSDPInvertS(DSDP dsdp)
Invert the S variables in each cone.
Definition: dsdpcops.c:307
int DSDPComputeXVariables(DSDP dsdp, double xmakermu, DSDPVec xmakery, DSDPVec xmakerdy, DSDPVec AX, double *tracexs)
Compute the X variables in each cone.
Definition: dsdpcops.c:654
int DSDPAddSchurRow(DSDP, int, DSDPVec)
Add a row to the Schur matrix.
Definition: dsdpcops.c:622
int DSDPComputeSS(DSDP dsdp, DSDPVec Y, DSDPDualFactorMatrix flag, DSDPTruth *ispsdefinite)
Compute the dual variables S in each cone.
Definition: dsdpcops.c:272
int DSDPComputeANorm2(DSDP dsdp, DSDPVec Anorm2)
Compute norm of A and C.
Definition: dsdpcops.c:246
int DSDPMonitorCones(DSDP dsdp, int tag)
This routine is called once per iteration.
Definition: dsdpcops.c:450
int DSDPSetSchurMatOps(DSDP, struct DSDPSchurMat_Ops *, void *)
Set the Schur complement matrix.
Definition: dsdpcops.c:602
int DSDPSetUpCones2(DSDP dsdp, DSDPVec yy0, DSDPSchurMat M)
Each cone should allocate its data structures .
Definition: dsdpcops.c:84
int DSDPGetConicDimension(DSDP dsdp, double *n)
Get the total dimension of the cones.
Definition: dsdpcops.c:401
int DSDPHessianMultiplyAdd(DSDP dsdp, DSDPVec v, DSDPVec vv)
Add the product of Schur matrix with v.
Definition: dsdpcops.c:188
int DSDPComputeG(DSDP dsdp, DSDPVec vt, DSDPVec vrhs1, DSDPVec vrhs2)
Compute the gradient of the barrier for each cone.
Definition: dsdpcops.c:215
int DSDPSchurMatReducePVec(DSDPSchurMat M, DSDPVec x)
Collect elements of the vector.
Definition: dsdpschurmat.c:307
int DSDPSchurMatSetR(DSDPSchurMat M, double rr)
Set up the data structure.
Definition: dsdpschurmat.c:338
int DSDPSchurMatZeroEntries(DSDPSchurMat M)
Zero all element in the matrix.
Definition: dsdpschurmat.c:97
int DSDPSchurMatAssemble(DSDPSchurMat M)
Final assembly of M.
Definition: dsdpschurmat.c:174
int DSDPSchurMatSetData(DSDPSchurMat *M, struct DSDPSchurMat_Ops *ops, void *data)
Set the Schur matrix with an opaque pointer and structure of function pointers.
Definition: dsdpschurmat.c:28
int DSDPSchurMatRowScaling(DSDPSchurMat M, DSDPVec D)
Identify which rows on on this processor.
Definition: dsdpschurmat.c:399
int DSDPSchurMatAddRow(DSDPSchurMat, int, double, DSDPVec)
Add elements to a row of the Schur matrix.
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
Schur complement matrix whose solution is the Newton direction.
Definition: dsdpschurmat.h:35
Internal structures for the DSDP solver.
Definition: dsdp.h:65