DSDP
sdpkcone.c
Go to the documentation of this file.
1#include "dsdpsdp.h"
2#include "dsdpcone_impl.h"
3#include "dsdpsys.h"
7static int SDPConeOperationsInitialize(struct DSDPCone_Ops*);
8
9static int KSDPConeSetup(void*,DSDPVec);
10static int KSDPConeSetup2(void*, DSDPVec, DSDPSchurMat);
11static int KSDPConeSize(void*,double*);
12static int KSDPConeSparsity(void*,int,int[],int[],int);
13static int KSDPConeComputeHessian(void*,double,DSDPSchurMat,DSDPVec,DSDPVec);
14static int KSDPConeComputeMaxStepLength(void*, DSDPVec, DSDPDualFactorMatrix, double *);
15static int KSDPConeComputeSS(void*, DSDPVec, DSDPDualFactorMatrix, DSDPTruth *);
16static int KSDPConeComputeLogSDeterminant(void *, double *, double*);
17static int KSDPConeComputeXX(void*, double, DSDPVec,DSDPVec,DSDPVec,double*);
18static int KSDPConeDestroy(void*);
19
20#undef __FUNCT__
21#define __FUNCT__ "KSDPConeComputeHessian"
22static int KSDPConeComputeHessian( void *K, double mu, DSDPSchurMat M, DSDPVec vrhs1, DSDPVec vrhs2){
23 int info;
24 SDPCone sdpcone=(SDPCone)K;
25 DSDPFunctionBegin;
26 info=SDPConeComputeHessian(sdpcone,mu,M,vrhs1,vrhs2);DSDPCHKERR(info);
27 DSDPFunctionReturn(0);
28}
29
30#undef __FUNCT__
31#define __FUNCT__ "KDPConeMultiply"
32static int KSDPConeMultiply( void *K, double mu, DSDPVec vrow, DSDPVec vin, DSDPVec vout){
33 int kk,info;
34 SDPCone sdpcone=(SDPCone)K;
35 SDPConeValid(sdpcone);
36 DSDPFunctionBegin;
37 for (kk=0; kk<sdpcone->nblocks; kk++){
38 info=SDPConeMultiply(sdpcone,kk,mu,vrow,vin,vout);DSDPCHKBLOCKERR(kk,info);
39 }
40 DSDPFunctionReturn(0);
41}
42
43#undef __FUNCT__
44#define __FUNCT__ "KDPConeRHS "
45static int KSDPConeRHS( void *K, double mu, DSDPVec vrow, DSDPVec vrhs1, DSDPVec vrhs2){
46 int kk,info;
47 SDPCone sdpcone=(SDPCone)K;
48 DSDPFunctionBegin;
49 SDPConeValid(sdpcone);
50 for (kk=0; kk<sdpcone->nblocks; kk++){
51 if (sdpcone->blk[kk].n<1) continue;
52 info=SDPConeComputeRHS(sdpcone,kk,mu,vrow,vrhs1,vrhs2); DSDPCHKBLOCKERR(kk,info);
53 }
54 DSDPFunctionReturn(0);
55}
56
57
58#undef __FUNCT__
59#define __FUNCT__ "KSDPConeSetup"
60static int KSDPConeSetup(void* K, DSDPVec y){
61 int info;
62 SDPCone sdpcone=(SDPCone)K;
63 DSDPFunctionBegin;
64 SDPConeValid(sdpcone);
65 info=SDPConeSetup(sdpcone,y);DSDPCHKERR(info);
66 DSDPFunctionReturn(0);
67}
68
69#undef __FUNCT__
70#define __FUNCT__ "KSDPConeSetup2"
71static int KSDPConeSetup2(void* K, DSDPVec y, DSDPSchurMat M){
72 int info;
73 SDPCone sdpcone=(SDPCone)K;
74 DSDPFunctionBegin;
75 info=SDPConeSetup2(sdpcone,y,M); DSDPCHKERR(info);
76 DSDPFunctionReturn(0);
77}
78
79
80#undef __FUNCT__
81#define __FUNCT__ "KSDPConeDestroy"
82static int KSDPConeDestroy(void* K){
83 int info;
84 SDPCone sdpcone=(SDPCone)K;
85 DSDPFunctionBegin;
86 SDPConeValid(sdpcone);
87 info=SDPConeDestroy(sdpcone);DSDPCHKERR(info);
88 DSDPFunctionReturn(0);
89}
90
91
92#undef __FUNCT__
93#define __FUNCT__ "KSDPConeSize"
94static int KSDPConeSize(void* K,double *n){
95 SDPCone sdpcone=(SDPCone)K;
96 DSDPFunctionBegin;
97 SDPConeValid(sdpcone);
98 *n=sdpcone->nn;
99 DSDPFunctionReturn(0);
100}
101
102#undef __FUNCT__
103#define __FUNCT__ "KSDPConeSparsity"
104static int KSDPConeSparsity(void *K,int row, int *tnnz, int rnnz[], int m){
105 SDPCone sdpcone=(SDPCone)K;
106 SDPblk *blk=sdpcone->blk;
107 int info,j,kk;
108 int nnzblocks=sdpcone->ATR.nnzblocks[row],*nzblocks=sdpcone->ATR.nzblocks[row];
109 DSDPFunctionBegin;
110 SDPConeValid(sdpcone);
111 for (j=0; j<nnzblocks; j++){
112 kk=nzblocks[j];
113 if (blk[kk].n<1) continue;
114 info=DSDPBlockDataMarkNonzeroMatrices(&blk[kk].ADATA,rnnz);DSDPCHKBLOCKERR(kk,info);
115 }
116 DSDPFunctionReturn(0);
117}
118
119#undef __FUNCT__
120#define __FUNCT__ "KSDPConeComputeSS"
121static int KSDPConeComputeSS(void *K, DSDPVec Y, DSDPDualFactorMatrix flag, DSDPTruth *ispsdefinite){
122 int kk,info;
123 SDPCone sdpcone=(SDPCone)K;
124 SDPblk *blk=sdpcone->blk;
125 DSDPTruth psdefinite;
126 DSDPDualMat SS;
127
128 DSDPFunctionBegin;
129 psdefinite = DSDP_TRUE;
130 for (kk=sdpcone->nblocks-1; kk>=0 && psdefinite == DSDP_TRUE; kk--){
131 if (blk[kk].n<1) continue;
132 if (flag==DUAL_FACTOR) SS=blk[kk].S;
133 else SS=blk[kk].SS;
134
135 switch (sdpcone->optype){
136 default:
137 info=SDPConeComputeSS(sdpcone,kk,Y,blk[kk].T);DSDPCHKBLOCKERR(kk,info);
138 info=DSDPDualMatSetArray(SS,blk[kk].T); DSDPCHKBLOCKERR(kk,info);
139 info=DSDPDualMatCholeskyFactor(SS,&psdefinite); DSDPCHKBLOCKERR(kk,info);
140 if (psdefinite == DSDP_FALSE){
141 if (flag==DUAL_FACTOR){
142 DSDPLogInfo(0,2,"Dual SDP Block %2.0f not PSD\n",kk);
143 } else {
144 DSDPLogInfo(0,2,"Primal SDP Block %2.0f not PSD\n",kk);
145 }
146 }
147 break;
148 }
149 }
150 *ispsdefinite=psdefinite;
151 if (flag==DUAL_FACTOR && psdefinite==DSDP_TRUE){
152 info=DSDPVecCopy(Y,sdpcone->YY);DSDPCHKERR(info);
153 }
154 DSDPFunctionReturn(0);
155}
156
157#undef __FUNCT__
158#define __FUNCT__ "KSDPConeInvertSS"
159static int KSDPConeInvertSS(void *K){
160 int kk,info;
161 SDPCone sdpcone=(SDPCone)K;
162 DSDPDualMat SS;
163
164 DSDPFunctionBegin;
165 SDPConeValid(sdpcone);
166 for (kk=0;kk<sdpcone->nblocks;kk++){
167 if (sdpcone->blk[kk].n<1) continue;
168 SS=sdpcone->blk[kk].S;
169 info=DSDPDualMatInvert(SS);DSDPCHKBLOCKERR(kk,info);
170 }
171 DSDPFunctionReturn(0);
172}
173
174
175
176#undef __FUNCT__
177#define __FUNCT__ "KSDPConeComputeMaxStepLength"
178static int KSDPConeComputeMaxStepLength(void *K, DSDPVec DY, DSDPDualFactorMatrix flag, double *maxsteplength){
179 int kk,info;
180 double smaxstep,maxmaxstep=1.0e20;
181 SDPCone sdpcone=(SDPCone)K;
182 DSDPDualMat SS;
183 SDPblk *blk=sdpcone->blk;
184 DSDPDSMat DS;
185 DSDPVMat T;
186
187 DSDPFunctionBegin;
188 SDPConeValid(sdpcone);
189 for (kk=0; kk<sdpcone->nblocks; kk++){
190 if (blk[kk].n<1) continue;
191 if (flag==DUAL_FACTOR) SS=blk[kk].S;
192 else SS=blk[kk].SS;
193 DS=blk[kk].DS; T=blk[kk].T;
194
195 info=SDPConeComputeSS(sdpcone,kk,DY,T);DSDPCHKBLOCKERR(kk,info);
196 info=DSDPDSMatSetArray(DS,T); DSDPCHKBLOCKERR(kk,info);
197
198 info=DSDPLanczosStepSize( &blk[kk].Lanczos,blk[kk].W,blk[kk].W2,SS,DS,&smaxstep );DSDPCHKBLOCKERR(kk,info);
199 DSDPLogInfo(0,12,"Block %d, PD %d, maxstepsize: %4.4e\n",kk,flag,smaxstep);
200 maxmaxstep=DSDPMin(smaxstep,maxmaxstep);
201 }
202 *maxsteplength=maxmaxstep;
203 DSDPFunctionReturn(0);
204}
205
206
207
208#undef __FUNCT__
209#define __FUNCT__ "KSDPConeAddANorm2"
210static int KSDPConeAddANorm2(void *K, DSDPVec ANorm2){
211 int kk,info;
212 SDPCone sdpcone=(SDPCone)K;
213 SDPblk *blk=sdpcone->blk;
214
215 DSDPFunctionBegin;
216 SDPConeValid(sdpcone);
217 for (kk=0; kk<sdpcone->nblocks; kk++){
218 if (blk[kk].n<1) continue;
219 info=DSDPBlockANorm2( &blk[kk].ADATA,ANorm2,blk[kk].n); DSDPCHKBLOCKERR(kk,info);
220 }
221 DSDPFunctionReturn(0);
222}
223
224
225
226#undef __FUNCT__
227#define __FUNCT__ "KSDPConeSetX"
228static int KSDPConeSetX(void *K, double mu, DSDPVec Y,DSDPVec DY){
229 SDPCone sdpcone=(SDPCone)K;
230 int info;
231 DSDPFunctionBegin;
232 SDPConeValid(sdpcone);
233 info=DSDPVecCopy(Y,sdpcone->YX);DSDPCHKERR(info);
234 info=DSDPVecCopy(DY,sdpcone->DYX);DSDPCHKERR(info);
235 sdpcone->xmakermu=mu;
236 DSDPFunctionReturn(0);
237}
238
239
240#undef __FUNCT__
241#define __FUNCT__ "KSDPConeComputeXX"
242static int KSDPConeComputeXX(void *K, double mu, DSDPVec Y,DSDPVec DY,DSDPVec AX,double* tracexs){
243
244 SDPCone sdpcone=(SDPCone)K;
245 int info,kk;
246 double xnorm,trxs,xtrace;
247 DSDPVMat X;
248
249 DSDPFunctionBegin;
250 SDPConeValid(sdpcone);
251 info=KSDPConeSetX(K,mu,Y,DY);DSDPCHKERR(info);
252 for (kk=0; kk<sdpcone->nblocks; kk++){
253 if (sdpcone->blk[kk].n<1) continue;
254 X=sdpcone->blk[kk].T;
255 info=SDPConeComputeX3(sdpcone,kk,mu,Y,DY,X);DSDPCHKBLOCKERR(kk,info);
256 info=SDPConeComputeXDot(sdpcone,kk,Y,X,AX,&xtrace,&xnorm,&trxs);DSDPCHKBLOCKERR(kk,info);
257 *tracexs+=trxs;
258 DSDPLogInfo(0,10,"SDP Block %d: norm(X): %4.2e, trace(X): %4.2e, trace(XS): %4.2e.\n",kk,xnorm,xtrace,trxs);
259 }
260 DSDPFunctionReturn(0);
261}
262
263
264#undef __FUNCT__
265#define __FUNCT__ "KSDPConeComputeLogSDeterminant"
266static int KSDPConeComputeLogSDeterminant(void *K, double *logdetobj, double *logdet){
267 int kk,info;
268 double dlogdet=0,dlogdet2=0,dd;
269 SDPCone sdpcone=(SDPCone)K;
270 SDPblk *blk=sdpcone->blk;
271
272 DSDPFunctionBegin;
273 SDPConeValid(sdpcone);
274 for (kk=0; kk<sdpcone->nblocks; kk++){
275 if (blk[kk].n<1) continue;
276 info=DSDPDualMatLogDeterminant(blk[kk].S,&dd);DSDPCHKBLOCKERR(kk,info);
277 dlogdet+=dd*blk[kk].gammamu;
278 dlogdet2+=dd*blk[kk].bmu;
279 }
280 *logdet=dlogdet;
281 *logdetobj=dlogdet2;
282 DSDPFunctionReturn(0);
283}
284
285
286#undef __FUNCT__
287#define __FUNCT__ "KSDPConeMonitor"
288int KSDPConeMonitor(void *K, int tag){
289 DSDPFunctionBegin;
290 DSDPFunctionReturn(0);
291}
292
293static struct DSDPCone_Ops kops;
294static const char *sdpconename ="SDP Cone";
295
296#undef __FUNCT__
297#define __FUNCT__ "SDPConeOperationsInitialize"
298static int SDPConeOperationsInitialize(struct DSDPCone_Ops* coneops){
299 int info;
300 if (coneops==NULL) return 0;
301 info=DSDPConeOpsInitialize(coneops); DSDPCHKERR(info);
302 coneops->conehessian=KSDPConeComputeHessian;
303 coneops->conerhs=KSDPConeRHS;
304 coneops->conesetup=KSDPConeSetup;
305 coneops->conesetup2=KSDPConeSetup2;
306 coneops->conedestroy=KSDPConeDestroy;
307 coneops->conecomputes=KSDPConeComputeSS;
308 coneops->coneinverts=KSDPConeInvertSS;
309 coneops->conesetxmaker=KSDPConeSetX;
310 coneops->conecomputex=KSDPConeComputeXX;
311 coneops->conemaxsteplength=KSDPConeComputeMaxStepLength;
312 coneops->conelogpotential=KSDPConeComputeLogSDeterminant;
313 coneops->conesize=KSDPConeSize;
314 coneops->conesparsity=KSDPConeSparsity;
315 coneops->conehmultiplyadd=KSDPConeMultiply;
316 coneops->coneanorm2=KSDPConeAddANorm2;
317 coneops->conemonitor=KSDPConeMonitor;
318 coneops->id=1;
319 coneops->name=sdpconename;
320 return 0;
321}
322
323#undef __FUNCT__
324#define __FUNCT__ "DSDPAddSDP"
331int DSDPAddSDP(DSDP dsdp,SDPCone sdpcone){
332 int info;
333 DSDPFunctionBegin;
334 SDPConeValid(sdpcone);
335 info=SDPConeOperationsInitialize(&kops); DSDPCHKERR(info);
336 info=DSDPAddCone(dsdp,&kops,(void*)sdpcone); DSDPCHKERR(info);
337 DSDPFunctionReturn(0);
338}
339
struct SDPCone_C * SDPCone
The SDPCone object points to blocks of data that specify semidefinite matrix inequalities.
Definition: dsdp5.h:26
DSDPDualFactorMatrix
DSDP requires two instances of the data structures S.
@ DUAL_FACTOR
DSDPTruth
Boolean variables.
@ DSDP_FALSE
@ DSDP_TRUE
int DSDPBlockDataMarkNonzeroMatrices(DSDPBlockData *ADATA, int *annz)
Mark which variable in block have a data matrix.
Definition: dsdpblock.c:254
int DSDPConeOpsInitialize(struct DSDPCone_Ops *dops)
Initialize the function pointers to 0.
Definition: dsdpcone.c:443
Implementations of a cone (SDP,LP,...) must provide a structure of function pointers.
int DSDPAddCone(DSDP, struct DSDPCone_Ops *, void *)
Apply DSDP to a conic structure.
Definition: dsdpcops.c:569
int DSDPDSMatSetArray(DSDPDSMat A, DSDPVMat T)
Set values into the matrix.
Definition: dsdpdsmat.c:130
int DSDPDualMatCholeskyFactor(DSDPDualMat S, DSDPTruth *psdefinite)
Factor the matrix.
Definition: dsdpdualmat.c:320
int DSDPDualMatLogDeterminant(DSDPDualMat S, double *logdet)
Free the matrix structure.
Definition: dsdpdualmat.c:122
int DSDPDualMatSetArray(DSDPDualMat S, DSDPVMat T)
Print the matrix.
Definition: dsdpdualmat.c:160
int DSDPDualMatInvert(DSDPDualMat S)
Invert the matrix.
Definition: dsdpdualmat.c:186
int DSDPLanczosStepSize(DSDPLanczosStepLength *, SDPConeVec, SDPConeVec, DSDPDualMat, DSDPDSMat, double *)
Compute distance to boundary.
Definition: dsdpstep.c:247
Internal SDPCone data structures and routines.
int SDPConeMultiply(SDPCone, int, double, DSDPVec, DSDPVec, DSDPVec)
Compute the gradient to the barrier term.
Definition: sdpcompute.c:182
int SDPConeComputeX3(SDPCone, int, double, DSDPVec, DSDPVec, DSDPVMat)
Compute the matrix X with the given information.
Definition: sdpcone.c:140
int SDPConeComputeRHS(SDPCone, int, double, DSDPVec, DSDPVec, DSDPVec)
Compute the gradient to the barrier term.
Definition: sdpcompute.c:125
int SDPConeComputeHessian(SDPCone, double, DSDPSchurMat, DSDPVec, DSDPVec)
Compute the Hessian to the barrier term.
Definition: sdpcompute.c:30
int SDPConeSetup2(SDPCone, DSDPVec, DSDPSchurMat)
Allocate data structure of the cone.
Definition: sdpconesetup.c:224
int SDPConeDestroy(SDPCone)
Free data structure of the cone.
Definition: sdpconesetup.c:350
int SDPConeSetup(SDPCone, DSDPVec)
Allocate data structure of the cone.
Definition: sdpconesetup.c:249
int SDPConeComputeSS(SDPCone, int, DSDPVec, DSDPVMat)
Sum the data matrices.
Definition: sdpcone.c:18
int SDPConeComputeXDot(SDPCone, int, DSDPVec, DSDPVMat, DSDPVec, double *, double *, double *)
Compute inner product of X with the Data, S, and norm of X.
Definition: sdpcone.c:111
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 DSDPAddSDP(DSDP dsdp, SDPCone sdpcone)
Pass a semidefinite cone to the solver.
Definition: sdpkcone.c:331
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
Schur complement matrix whose solution is the Newton direction.
Definition: dsdpschurmat.h:35
Dense symmetric matrix for one block in the semidefinite cone.
Definition: dsdpxmat.h:17
Internal structures for the DSDP solver.
Definition: dsdp.h:65
Internal structure for semidefinite cone.
Definition: dsdpsdp.h:80
Internal structure for block of semidefinite cone.
Definition: dsdpsdp.h:52