DSDP
sdpconesetup.c
Go to the documentation of this file.
1#include "dsdpsdp.h"
2#include "dsdpsys.h"
8#undef __FUNCT__
9#define __FUNCT__ "DSDPDataTransposeInitialize"
16 DSDPFunctionBegin;
17 ATranspose->nnzblocks=0;
18 ATranspose->nzblocks=0;
19 ATranspose->idA=0;
20 ATranspose->idAP=0;
21 ATranspose->ttnzmat=0;
22 ATranspose->nnzblocks=0;
23 DSDPFunctionReturn(0);
24}
25
26#undef __FUNCT__
27#define __FUNCT__ "DSDPDataTransposeSetup"
36int DSDPDataTransposeSetup(DSDPDataTranspose *ATranspose, SDPblk *blk, int nblocks, int m){
37
38 int i,ii,kk,vvar,info;
39 int nnzmats,tnzmats=0;
40 DSDPFunctionBegin;
41
42 info=DSDPDataTransposeTakeDown(ATranspose);DSDPCHKERR(info);
43 /* Determine sparsity pattern of SDP Data Matrices */
44
45 DSDPCALLOC2(&ATranspose->nnzblocks,int,(m),&info);DSDPCHKERR(info);
46 DSDPCALLOC2(&ATranspose->nzblocks,int*,(m),&info);DSDPCHKERR(info);
47 DSDPCALLOC2(&ATranspose->idA,int*,(m),&info);DSDPCHKERR(info);
48 ATranspose->m=m;
49 for (i=0;i<m;i++){ ATranspose->nnzblocks[i]=0; }
50 for (kk=0; kk<nblocks; kk++){
51 info=DSDPBlockDataMarkNonzeroMatrices(&blk[kk].ADATA,ATranspose->nnzblocks);DSDPCHKERR(info);
52 }
53 for (tnzmats=0,i=0;i<m;i++){ tnzmats += ATranspose->nnzblocks[i];}
54
55 DSDPCALLOC2(&ATranspose->ttnzmat,int,tnzmats,&info);DSDPCHKERR(info);
56 ATranspose->nzblocks[0]=ATranspose->ttnzmat;
57 for (i=1;i<m;i++){
58 ATranspose->nzblocks[i]=ATranspose->nzblocks[i-1]+ATranspose->nnzblocks[i-1];
59 }
60
61 DSDPCALLOC2(&ATranspose->idAP,int,tnzmats,&info);DSDPCHKERR(info);
62 ATranspose->idA[0]=ATranspose->idAP;
63 for (i=1;i<m;i++){
64 ATranspose->idA[i]=ATranspose->idA[i-1]+ATranspose->nnzblocks[i-1];
65 }
66
67 for (i=0;i<m;i++){ATranspose->nnzblocks[i]=0;}
68 for (kk=0; kk<nblocks; kk++){
69 info=DSDPBlockCountNonzeroMatrices(&blk[kk].ADATA,&nnzmats);DSDPCHKERR(info);
70 for (i=0;i<nnzmats;i++){
71 info=DSDPBlockGetMatrix(&blk[kk].ADATA,i,&ii,0,0);DSDPCHKERR(info);
72 vvar=ATranspose->nnzblocks[ii];
73 ATranspose->nzblocks[ii][vvar]=kk;
74 ATranspose->idA[ii][vvar]=i;
75 ATranspose->nnzblocks[ii]++;
76 }
77 }
78
79 DSDPFunctionReturn(0);
80}
81
82#undef __FUNCT__
83#define __FUNCT__ "DSDPDataTransposeTakeDown"
90 int info;
91 DSDPFunctionBegin;
92 DSDPFREE(&ATranspose->ttnzmat,&info);DSDPCHKERR(info);
93 DSDPFREE(&ATranspose->idAP,&info);DSDPCHKERR(info);
94 DSDPFREE(&ATranspose->nzblocks,&info);DSDPCHKERR(info);
95 DSDPFREE(&ATranspose->nnzblocks,&info);DSDPCHKERR(info);
96 DSDPFREE(&ATranspose->idA,&info);DSDPCHKERR(info);
97 info=DSDPDataTransposeInitialize(ATranspose);DSDPCHKERR(info);
98 DSDPFunctionReturn(0);
99}
100
101#undef __FUNCT__
102#define __FUNCT__ "DSDPCreateSDPCone"
113int DSDPCreateSDPCone(DSDP dsdp, int blocks, SDPCone* dspcone){
114 int i,info;
115 SDPCone sdpcone;
116
117 DSDPFunctionBegin;
118 DSDPCALLOC1(&sdpcone,struct SDPCone_C,&info);DSDPCHKERR(info);
119 *dspcone=sdpcone;
120 sdpcone->keyid=SDPCONEKEY;
121 info=DSDPAddSDP(dsdp,sdpcone);DSDPCHKERR(info);
122
123 info=DSDPGetNumberOfVariables(dsdp,&sdpcone->m);DSDPCHKERR(info);
124 DSDPCALLOC2(&sdpcone->blk,SDPblk,blocks,&info); DSDPCHKERR(info);
125 for (i=0;i<blocks; i++){
126 info=DSDPBlockInitialize(&sdpcone->blk[i]); DSDPCHKBLOCKERR(i,info);
127 }
128
129 sdpcone->nblocks=blocks;
130 sdpcone->optype=3;
131 info=DSDPUseDefaultDualMatrix(sdpcone); DSDPCHKERR(info);
132
133 sdpcone->nn=0;
134 sdpcone->dsdp=dsdp;
135 info=DSDPDataTransposeInitialize(&sdpcone->ATR); DSDPCHKERR(info);
136 info=DSDPBlockEventZero();DSDPCHKERR(info);
137 info=DSDPDualMatEventZero();DSDPCHKERR(info);
138 info=DSDPVMatEventZero();DSDPCHKERR(info);
139 DSDPFunctionReturn(0);
140}
141
142
144
145#undef __FUNCT__
146#define __FUNCT__ "DSDPBlockSetup"
154int DSDPBlockSetup(SDPblk *blk, int blockj, DSDPVec WY){
155 int n,info,trank,flag;
156 DSDPFunctionBegin;
157 /*
158 info=DSDPBlockTakeDown(blk);DSDPCHKERR(info);
159 */
160 n=blk->n;
161 info=DSDPVMatExist(blk->T,&flag);DSDPCHKERR(info);
162 if (flag==0){
163 info=DSDPMakeVMat(blk->format,n,&blk->T);DSDPCHKERR(info);
164 }
165
166 info = DSDPIndexCreate(blk->n,&blk->IS);DSDPCHKERR(info);
167 info = SDPConeVecCreate(blk->n,&blk->W);DSDPCHKERR(info);
168 info = SDPConeVecDuplicate(blk->W,&blk->W2);DSDPCHKERR(info);
169
170 /* Build Lanczos structures */
171 info=DSDPSetMaximumLanczosIterations(&blk->Lanczos,20); DSDPCHKERR(info);
172 if (n>30){info=DSDPSetMaximumLanczosIterations(&blk->Lanczos,20); DSDPCHKERR(info);}
173 if (n>300){info=DSDPSetMaximumLanczosIterations(&blk->Lanczos,40); DSDPCHKERR(info);}
174 if (n>1000){info=DSDPSetMaximumLanczosIterations(&blk->Lanczos,50); DSDPCHKERR(info);}
175 if (1){
176 info=DSDPFastLanczosSetup(&blk->Lanczos,blk->W);DSDPCHKERR(info);
177 DSDPLogInfo(0,19,"SDP Block %d using Fast Lanczos\n",blockj);
178 } else {
179 info=DSDPRobustLanczosSetup(&blk->Lanczos,blk->W);DSDPCHKERR(info);
180 DSDPLogInfo(0,19,"SDP Block %d using Full Lanczos\n",blockj);
181 }
182
183 /* Eigenvalues and Eigenvectors */
184 info=DSDPBlockFactorData(&blk->ADATA,blk->T,blk->W);DSDPCHKERR(info);
185 info=DSDPBlockDataRank(&blk->ADATA,&trank,n);DSDPCHKERR(info);
186
187 info=DSDPCreateS(&blk->ADATA,blk->format,trank,WY,blk->T,blk->W,blk->W2,&blk->S,&blk->SS,&blk->DS,0);DSDPCHKERR(info);
188
189 DSDPFunctionReturn(0);
190}
191
192#undef __FUNCT__
193#define __FUNCT__ "SDPConeBlockNNZ"
194int SDPConeBlockNNZ(SDPblk *blk,int m){
195 int i,ii,n,info,nnz,nnzmats,tnnzmats,tnnz=0;
196 double scl;
197 DSDPDataMat AA;
198 DSDPFunctionBegin;
199 nnzmats=blk->ADATA.nnzmats;tnnzmats=nnzmats;
200 n=blk->n;
201
202 for (i=0;i<nnzmats;i++){
203 info=DSDPBlockGetMatrix(&blk->ADATA,i,&ii,&scl,&AA);DSDPCHKERR(info);
204 if (ii==0){tnnzmats--; continue;}
205 if (ii==m-1){continue;}
206 info = DSDPDataMatCountNonzeros(AA,&nnz,n); DSDPCHKERR(info);
207 tnnz+= (nnz*(tnnzmats-i));
208 }
209 if (tnnzmats>1){ tnnz=tnnz/((tnnzmats)*(tnnzmats+1)/2); }
210 if (tnnz<1) tnnz = 1;
211 blk->nnz=tnnz;
212 DSDPFunctionReturn(0);
213}
214
215#undef __FUNCT__
216#define __FUNCT__ "SDPConeSetup2"
225 int kk,n,m,info;
226 double nn=0;
227 SDPblk *blk;
228 DSDPFunctionBegin;
229 info=DSDPVecGetSize(yy0,&m);DSDPCHKERR(info);
230 for (kk=0; kk<sdpcone->nblocks; kk++){
231 blk=&sdpcone->blk[kk];
232 n=blk->n;
233 info=SDPConeBlockNNZ(blk,m);DSDPCHKERR(info);
234 info=DSDPBlockSetup(blk,kk,sdpcone->Work);DSDPCHKERR(info);
235 nn+=n*blk->gammamu;
236 }
237 sdpcone->nn=(int)nn;
238 DSDPFunctionReturn(0);
239}
240
241#undef __FUNCT__
242#define __FUNCT__ "SDPConeSetup"
249int SDPConeSetup(SDPCone sdpcone, DSDPVec yy0){
250 int kk,n,m,info;
251 DSDPFunctionBegin;
252
253 info = DSDPVecGetSize(yy0,&m);DSDPCHKERR(info);
254 if (m!=sdpcone->m+2){DSDPSETERR(8,"CHECK DIMENSION\n");}
255 info = DSDPVecDuplicate(yy0,&sdpcone->Work);DSDPCHKERR(info);
256 info = DSDPVecDuplicate(yy0,&sdpcone->Work2);DSDPCHKERR(info);
257 info = DSDPVecDuplicate(yy0,&sdpcone->YY);DSDPCHKERR(info);
258 info = DSDPVecDuplicate(yy0,&sdpcone->YX);DSDPCHKERR(info);
259 info = DSDPVecDuplicate(yy0,&sdpcone->DYX);DSDPCHKERR(info);
260 for (kk=0; kk<sdpcone->nblocks; kk++){
261 n=sdpcone->blk[kk].n;
262 info=SDPConeSetRIdentity(sdpcone,kk,n,1.0);DSDPCHKERR(info);
263 }
264
265 info=DSDPDataTransposeSetup(&sdpcone->ATR,sdpcone->blk,sdpcone->nblocks,m); DSDPCHKERR(info);
266 info=DSDPBlockEventInitialize();DSDPCHKERR(info);
267 info=DSDPDualMatEventInitialize();DSDPCHKERR(info);
268 info=DSDPVMatEventInitialize();DSDPCHKERR(info);
269 DSDPFunctionReturn(0);
270}
271
272#undef __FUNCT__
273#define __FUNCT__ "DSDPBlockInitialize"
280 int info;
281 DSDPFunctionBegin;
282 blk->n=0;
283 blk->format='N';
284 blk->gammamu=1.0;
285 blk->bmu=0.0;
286 blk->nnz=10000000;
287
288 info = DSDPDualMatInitialize(&blk->S); DSDPCHKERR(info);
289 info = DSDPDualMatInitialize(&blk->SS); DSDPCHKERR(info);
290 info = DSDPDSMatInitialize(&blk->DS); DSDPCHKERR(info);
291 info = DSDPVMatInitialize(&blk->T); DSDPCHKERR(info);
292 info = DSDPLanczosInitialize(&blk->Lanczos); DSDPCHKERR(info);
293 info = DSDPBlockDataInitialize(&blk->ADATA); DSDPCHKERR(info);
294 info = DSDPIndexInitialize(&blk->IS); DSDPCHKERR(info);
295 DSDPFunctionReturn(0);
296}
297
298#undef __FUNCT__
299#define __FUNCT__ "DSDPBlockTakeDown"
306 int info;
307 DSDPFunctionBegin;
308 if (!blk){DSDPFunctionReturn(0);}
309 info=DSDPBlockTakeDownData(&blk->ADATA);DSDPCHKERR(info);
310 info=SDPConeVecDestroy(&blk->W);DSDPCHKERR(info);
311 info=SDPConeVecDestroy(&blk->W2);DSDPCHKERR(info);
312 info=DSDPIndexDestroy(&blk->IS);DSDPCHKERR(info);
313 info=DSDPLanczosDestroy(&blk->Lanczos);DSDPCHKERR(info);
314 info=DSDPDualMatDestroy(&blk->SS);DSDPCHKERR(info);
315 info=DSDPDualMatDestroy(&blk->S);DSDPCHKERR(info);
316 info=DSDPDSMatDestroy(&blk->DS);DSDPCHKERR(info);
317 info=DSDPVMatDestroy(&blk->T);DSDPCHKERR(info);
318 DSDPFunctionReturn(0);
319}
320
321#undef __FUNCT__
322#define __FUNCT__ "DSDPConeTakeDown"
329 int blockj,info;
330 DSDPFunctionBegin;
331 for (blockj=0; blockj<sdpcone->nblocks; blockj++){
332 info=DSDPBlockTakeDown(&sdpcone->blk[blockj]);DSDPCHKERR(info);
333 }
334 info=DSDPVecDestroy(&sdpcone->Work);DSDPCHKERR(info);
335 info=DSDPVecDestroy(&sdpcone->Work2);DSDPCHKERR(info);
336 info=DSDPVecDestroy(&sdpcone->YY);DSDPCHKERR(info);
337 info=DSDPVecDestroy(&sdpcone->YX);DSDPCHKERR(info);
338 info=DSDPVecDestroy(&sdpcone->DYX);DSDPCHKERR(info);
339 info=DSDPDataTransposeTakeDown(&sdpcone->ATR);DSDPCHKERR(info);
340 DSDPFunctionReturn(0);
341}
342
343#undef __FUNCT__
344#define __FUNCT__ "SDPConeDestroy"
351 int blockj,info;
352 DSDPFunctionBegin;
353 info=DSDPConeTakeDown(sdpcone);DSDPCHKERR(info);
354 for (blockj=0; blockj<sdpcone->nblocks; blockj++){
355 info=DSDPBlockDataDestroy(&sdpcone->blk[blockj].ADATA);DSDPCHKERR(info);
356 }
357 DSDPFREE(&sdpcone->blk,&info);DSDPCHKERR(info);
358 DSDPFREE(&sdpcone,&info);DSDPCHKERR(info);
359 info=DSDPBlockEventZero();DSDPCHKERR(info);
360 info=DSDPDualMatEventZero();DSDPCHKERR(info);
361 info=DSDPVMatEventZero();DSDPCHKERR(info);
362 DSDPFunctionReturn(0);
363}
364
int SDPConeSetRIdentity(SDPCone sdpcone, int blockj, int n, double rr)
Add identify matrix to dual matrix.
int DSDPBlockDataDestroy(DSDPBlockData *ADATA)
Free the data matrices.
Definition: dsdpblock.c:195
int DSDPBlockGetMatrix(DSDPBlockData *ADATA, int id, int *vari, double *scl, DSDPDataMat *A)
Get a data matrix from a block of data.
Definition: dsdpblock.c:307
int DSDPBlockDataInitialize(DSDPBlockData *ADATA)
Set pointers to null.
Definition: dsdpblock.c:163
int DSDPBlockDataMarkNonzeroMatrices(DSDPBlockData *ADATA, int *annz)
Mark which variable in block have a data matrix.
Definition: dsdpblock.c:254
int DSDPBlockCountNonzeroMatrices(DSDPBlockData *ADATA, int *nzmats)
Count how many data matrices are in a block of data.
Definition: dsdpblock.c:272
int DSDPBlockFactorData(DSDPBlockData *ADATA, DSDPVMat X, SDPConeVec W)
Factor the data matrices.
Definition: dsdpblock.c:113
int DSDPBlockTakeDownData(DSDPBlockData *ADATA)
Free structures in block of data.
Definition: dsdpblock.c:182
int DSDPDataMatCountNonzeros(DSDPDataMat A, int *nnz, int n)
Compute the square of the Frobenius norm.
Definition: dsdpdatamat.c:152
int DSDPDSMatDestroy(DSDPDSMat *A)
Free the data structure.
Definition: dsdpdsmat.c:70
int DSDPDSMatInitialize(DSDPDSMat *B)
Set pointers to null.
Definition: dsdpdsmat.c:254
int DSDPDualMatInitialize(DSDPDualMat *S)
Set pointers to null.
Definition: dsdpdualmat.c:471
int DSDPDualMatDestroy(DSDPDualMat *S)
Free the matrix structure.
Definition: dsdpdualmat.c:65
int DSDPFastLanczosSetup(DSDPLanczosStepLength *, SDPConeVec)
Use Lanczos procedure. Assume off tridiagonal entries are zero.
Definition: dsdpstep.c:133
int DSDPLanczosDestroy(DSDPLanczosStepLength *)
Free data structure.
Definition: dsdpstep.c:191
int DSDPLanczosInitialize(DSDPLanczosStepLength *)
Initialize Lanczos structure.
Definition: dsdpstep.c:92
int DSDPSetMaximumLanczosIterations(DSDPLanczosStepLength *LZ, int)
Set parameter.
Definition: dsdpstep.c:119
int DSDPRobustLanczosSetup(DSDPLanczosStepLength *, SDPConeVec)
Use slowerer but more robust method.
Definition: dsdpstep.c:163
Internal SDPCone data structures and routines.
int DSDPAddSDP(DSDP, SDPCone)
Pass a semidefinite cone to the solver.
Definition: sdpkcone.c:331
int DSDPMakeVMat(char, int, DSDPVMat *)
Allocate V matrix.
Definition: sdpsss.c:351
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 DSDPVMatInitialize(DSDPVMat *B)
Set pointers to null.
Definition: dsdpxmat.c:424
int DSDPVMatExist(DSDPVMat X, int *flag)
Answer whether the array has been allocated or not.
Definition: dsdpxmat.c:440
int DSDPVMatDestroy(DSDPVMat *X)
Deallocate matrix.
Definition: dsdpxmat.c:86
int DSDPGetNumberOfVariables(DSDP dsdp, int *m)
Copy the number of variables y.
Definition: dsdpsetdata.c:707
int DSDPConeTakeDown(SDPCone sdpcone)
Free data structure of the cone.
Definition: sdpconesetup.c:328
int DSDPBlockTakeDown(SDPblk *blk)
Free data structures in one block of the cone.
Definition: sdpconesetup.c:305
int DSDPBlockSetup(SDPblk *blk, int blockj, DSDPVec WY)
Allocate data structures of one block the cone.
Definition: sdpconesetup.c:154
int DSDPCreateS(DSDPBlockData *, char, int, DSDPVec, DSDPVMat, SDPConeVec, SDPConeVec, DSDPDualMat *, DSDPDualMat *, DSDPDSMat *, void *)
Create S1, S2, and DS.
Definition: sdpsss.c:314
int DSDPBlockInitialize(SDPblk *blk)
Initialize data structures in one block of the cone.
Definition: sdpconesetup.c:279
int SDPConeSetup(SDPCone sdpcone, DSDPVec yy0)
Allocate data structure of the cone.
Definition: sdpconesetup.c:249
int SDPConeDestroy(SDPCone sdpcone)
Free data structure of the cone.
Definition: sdpconesetup.c:350
int SDPConeSetup2(SDPCone sdpcone, DSDPVec yy0, DSDPSchurMat M)
Allocate data structure of the cone.
Definition: sdpconesetup.c:224
int DSDPDataTransposeSetup(DSDPDataTranspose *ATranspose, SDPblk *blk, int nblocks, int m)
Set up transpose structure for data.
Definition: sdpconesetup.c:36
int DSDPDataTransposeTakeDown(DSDPDataTranspose *ATranspose)
Free transpose structure for data.
Definition: sdpconesetup.c:89
int DSDPDataTransposeInitialize(DSDPDataTranspose *ATranspose)
Initialize transpose structure for data.
Definition: sdpconesetup.c:15
int SDPConeVecDuplicate(SDPConeVec V1, SDPConeVec *V2)
Allocate another vector with the same structure as the first.
Definition: sdpconevec.c:195
int DSDPIndexDestroy(DSDPIndex *IS)
Deallocate memory.
Definition: sdpconevec.c:264
int DSDPIndexCreate(int n, DSDPIndex *IS)
Allocate array for indices.
Definition: sdpconevec.c:248
int DSDPIndexInitialize(DSDPIndex *IS)
Set structure pointers to 0.
Definition: sdpconevec.c:234
Internal structure for data in one block of semidefintie.
Definition: dsdpsdp.h:39
Symmetric Delta S matrix for one block in the semidefinite cone.
Definition: dsdpdsmat.h:23
Symmetric data matrix for one block in the semidefinite cone.
Definition: dsdpdatamat.h:15
Internal structure for transpose of data.
Definition: dsdpsdp.h:25
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
Vector whose length corresponds to dimension of a block in a cone.
Definition: sdpconevec.h:13
Internal structure for semidefinite cone.
Definition: dsdpsdp.h:80
Internal structure for block of semidefinite cone.
Definition: dsdpsdp.h:52