DSDP
dsdpblock.c
Go to the documentation of this file.
1 #include "dsdpsdp.h"
2 #include "dsdpsys.h"
3 
4 static int sdpvecvecevent=0,sdpdotevent=0;
10 #undef __FUNCT__
11 #define __FUNCT__ "DSDPBlockASum"
12 
20 int DSDPBlockASum(DSDPBlockData *ADATA, double aa, DSDPVec Yk, DSDPVMat XX){
21 
22  double *xx,ytmp,scl=ADATA->scl;
23  int ii,vari,n,nn,info;
24 
25  DSDPFunctionBegin;
26  info=DSDPVMatGetSize(XX, &n); DSDPCHKERR(info);
27  info=DSDPVMatGetArray(XX, &xx, &nn); DSDPCHKERR(info);
28  for (ii=0;ii<ADATA->nnzmats;ii++){
29  vari=ADATA->nzmat[ii];
30  info=DSDPVecGetElement(Yk,vari,&ytmp);DSDPCHKVARERR(vari,info);
31  if (ytmp==0) continue;
32  info = DSDPDataMatAddMultiple(ADATA->A[ii], -aa*scl*ytmp, xx,nn,n); DSDPCHKVARERR(vari,info);
33  }
34  info=DSDPVMatRestoreArray(XX, &xx, &nn); DSDPCHKERR(info);
35  DSDPFunctionReturn(0);
36 }
37 
38 #undef __FUNCT__
39 #define __FUNCT__ "DSDPBlockADot"
40 
49 int DSDPBlockADot(DSDPBlockData *ADATA, double aa, DSDPVec Alpha, DSDPVMat X, DSDPVec AX){
50 
51  int ii,vari,n,nn,info;
52  double *x,sum=0,aalpha=0,scl=ADATA->scl;
53 
54  DSDPFunctionBegin;
55  DSDPEventLogBegin(sdpdotevent);
56  info=DSDPVMatScaleDiagonal(X,0.5); DSDPCHKERR(info);
57  info=DSDPVMatGetSize(X, &n); DSDPCHKERR(info);
58  info=DSDPVMatGetArray(X, &x, &nn); DSDPCHKERR(info);
59  for (ii=0;ii<ADATA->nnzmats; ii++){ /* Matrix Entries */
60  vari=ADATA->nzmat[ii];
61  info=DSDPVecGetElement(Alpha,vari,&aalpha);DSDPCHKVARERR(vari,info);
62  if (aalpha==0.0) continue;
63  info=DSDPDataMatDot(ADATA->A[ii],x,nn,n,&sum);DSDPCHKVARERR(vari,info);
64  info=DSDPVecAddElement(AX,vari,aa*aalpha*sum*scl);DSDPCHKVARERR(vari,info);
65  }
66  info=DSDPVMatRestoreArray(X, &x, &nn); DSDPCHKERR(info);
67  info=DSDPVMatScaleDiagonal(X,2.0); DSDPCHKERR(info);
68  DSDPEventLogEnd(sdpdotevent);
69  DSDPFunctionReturn(0);
70 }
71 
72 #undef __FUNCT__
73 #define __FUNCT__ "DSDPBlockvAv"
74 
84 int DSDPBlockvAv(DSDPBlockData *ADATA, double aa, DSDPVec Alpha, SDPConeVec V, DSDPVec VAV){
85 
86  int ii,vari,info;
87  double sum=0,aalpha=0,scl=ADATA->scl;
88 
89  DSDPFunctionBegin;
90  DSDPEventLogBegin(sdpvecvecevent);
91  if (aa==0){DSDPFunctionReturn(0);}
92  for (ii=0;ii<ADATA->nnzmats; ii++){ /* Matrix Entries */
93  vari=ADATA->nzmat[ii];
94  info=DSDPVecGetElement(Alpha,vari,&aalpha);DSDPCHKVARERR(vari,info);
95  if (aalpha==0.0) continue;
96  info=DSDPDataMatVecVec(ADATA->A[ii],V,&sum);DSDPCHKVARERR(vari,info);
97  info=DSDPVecAddElement(VAV,vari,aa*aalpha*sum*scl);DSDPCHKVARERR(vari,info);
98  }
99  DSDPEventLogEnd(sdpvecvecevent);
100  DSDPFunctionReturn(0);
101 }
102 
103 #undef __FUNCT__
104 #define __FUNCT__ "DSDPBlockFactorData"
105 
114  DSDPVMat X, SDPConeVec W){
115 
116  int ii,vari,n,nn,info,*iwork3n,i13,n26;
117  double *x,*dwork3n;
118  DSDPDataMat AA;
119 
120  DSDPFunctionBegin;
121  info=DSDPVMatGetSize(X, &n); DSDPCHKERR(info);
122  i13=13*n+1;n26=26*n+1;
123  DSDPCALLOC2(&dwork3n,double,n26,&info);DSDPCHKERR(info);
124  DSDPCALLOC2(&iwork3n,int,i13,&info);DSDPCHKERR(info);
125  info=DSDPVMatGetArray(X, &x, &nn); DSDPCHKERR(info);
126  for (ii=0;ii<ADATA->nnzmats; ii++){ /* Matrix Entries */
127  info=DSDPBlockGetMatrix(ADATA,ii,&vari,0,&AA);DSDPCHKVARERR(vari,info);
128  DSDPLogInfo(0,39,"SDP Data Mat Setup: %d\n",vari);
129  if (vari==0) continue;
130  info=DSDPDataMatFactor(AA,W,x,nn,dwork3n,n26,iwork3n,i13); DSDPCHKVARERR(vari,info);
131  }
132  info=DSDPVMatRestoreArray(X, &x, &nn); DSDPCHKERR(info);
133  DSDPFREE(&dwork3n,&info);DSDPCHKERR(info);
134  DSDPFREE(&iwork3n,&info);DSDPCHKERR(info);
135  DSDPFunctionReturn(0);
136 }
137 
138 
139 #undef __FUNCT__
140 #define __FUNCT__ "DSDPBlockEventZero"
141 int DSDPBlockEventZero(void){
142  DSDPFunctionBegin;
143  sdpvecvecevent=0;sdpdotevent=0;
144  DSDPFunctionReturn(0);
145 }
146 
147 #undef __FUNCT__
148 #define __FUNCT__ "DSDPBlockEventInitialize"
149 int DSDPBlockEventInitialize(void){
150  DSDPFunctionBegin;
151  if (sdpvecvecevent==0){DSDPEventLogRegister("SDP VecMatVec",&sdpvecvecevent);}
152  if (sdpdotevent==0){DSDPEventLogRegister("SDP Dot",&sdpdotevent);}
153  DSDPFunctionReturn(0);
154 }
155 
156 #undef __FUNCT__
157 #define __FUNCT__ "DSDPBlockDataInitialize"
158 
164  DSDPFunctionBegin;
165  ADATA->nnzmats=0;
166  ADATA->maxnnzmats=0;
167  ADATA->nzmat=0;
168  ADATA->A=0;
169  ADATA->r=1.0;
170  ADATA->scl=1.0;
171  /* ADATA->n=0; */
172  DSDPFunctionReturn(0);
173 }
174 
175 #undef __FUNCT__
176 #define __FUNCT__ "DSDPBlockTakeDownData"
177 
183  DSDPFunctionBegin;
184  DSDPFunctionReturn(0);
185 }
186 
187 
188 #undef __FUNCT__
189 #define __FUNCT__ "DSDPBlockDataDestroy"
190 
196  int ii,vari,info;
197  DSDPFunctionBegin;
198  if (!ADATA){DSDPFunctionReturn(0);}
199  DSDPLogInfo(0,18,"Destroying All Existing Data Matrices \n");
200  for (ii=0; ii<ADATA->nnzmats; ii++){
201  vari=ADATA->nzmat[ii];
202  info = DSDPDataMatDestroy(&ADATA->A[ii]);DSDPCHKVARERR(vari,info);
203  ADATA->nzmat[ii]=0;
204  }
205  ADATA->nnzmats=0;
206  info=DSDPBlockTakeDownData(ADATA);DSDPCHKERR(info);
207  DSDPFREE(&ADATA->nzmat,&info);DSDPCHKERR(info);
208  DSDPFREE(&ADATA->A,&info);DSDPCHKERR(info);
209  info=DSDPBlockDataInitialize(ADATA);DSDPCHKERR(info);
210  DSDPFunctionReturn(0);
211 }
212 
213 #undef __FUNCT__
214 #define __FUNCT__ "DSDPBlockDataAllocate"
215 
222  int j,info,*nzmat;
223  DSDPDataMat *A;
224  DSDPFunctionBegin;
225  if (!ADATA){DSDPFunctionReturn(0);}
226  if (nnz<=ADATA->maxnnzmats){DSDPFunctionReturn(0);}
227  DSDPLogInfo(0,18,"REALLOCATING SPACE FOR %d SDP BLOCK MATRICES! Previously allocated: %d \n",nnz,ADATA->maxnnzmats);
228  DSDPCALLOC2(&A,struct DSDPDataMat_C,nnz,&info);DSDPCHKERR(info);
229  DSDPCALLOC2(&nzmat,int,nnz,&info);DSDPCHKERR(info);
230  for (j=0;j<nnz;j++){nzmat[j]=0;}
231  for (j=0;j<nnz;j++){info = DSDPDataMatInitialize(&A[j]);DSDPCHKERR(info);}
232  if (ADATA->maxnnzmats>0){
233  for (j=0;j<ADATA->nnzmats;j++){nzmat[j]=ADATA->nzmat[j];}
234  for (j=0;j<ADATA->nnzmats;j++){A[j]=ADATA->A[j];}
235  DSDPFREE(&ADATA->A,&info);DSDPCHKERR(info);
236  DSDPFREE(&ADATA->nzmat,&info);DSDPCHKERR(info);
237  } else {
238  ADATA->nnzmats=0;
239  }
240  ADATA->maxnnzmats=nnz;
241  ADATA->nzmat=nzmat;
242  ADATA->A=A;
243  DSDPFunctionReturn(0);
244 }
245 
246 #undef __FUNCT__
247 #define __FUNCT__ "DSDPBlockDataMarkNonzeroMatrices"
248 
255  int i;
256  DSDPFunctionBegin;
257  for (i=0; i<ADATA->nnzmats; i++){
258  annz[ADATA->nzmat[i]]++;
259  }
260  DSDPFunctionReturn(0);
261 }
262 
263 #undef __FUNCT__
264 #define __FUNCT__ "DSDPBlockCountNonzerosMatrices"
265 
273  DSDPFunctionBegin;
274  *nzmats=ADATA->nnzmats;
275  DSDPFunctionReturn(0);
276 }
277 
278 #undef __FUNCT__
279 #define __FUNCT__ "DSDPBlockDataRank"
280 int DSDPBlockDataRank(DSDPBlockData *ADATA, int *trank, int n){
281  int ii,vari,info,ri,r2=0;
282  DSDPDataMat AA;
283 
284  DSDPFunctionBegin;
285  for (ii=0;ii<ADATA->nnzmats;ii++){
286  info=DSDPBlockGetMatrix(ADATA,ii,&vari,0,&AA);DSDPCHKVARERR(vari,info);
287  if (vari==0) continue;
288  info=DSDPDataMatGetRank(AA,&ri,n); DSDPCHKVARERR(vari,info);
289  r2+=ri;
290  }
291  *trank=r2;
292  DSDPFunctionReturn(0);
293 }
294 
295 #undef __FUNCT__
296 #define __FUNCT__ "DSDPBlockGetMatrix"
297 
307 int DSDPBlockGetMatrix(DSDPBlockData *ADATA,int id, int *vari, double *scl, DSDPDataMat *A){
308  DSDPFunctionBegin;
309  if (id>=0 && id < ADATA->nnzmats){
310  if (vari) *vari=ADATA->nzmat[id];
311  if (scl) *scl=ADATA->scl;
312  if (A) *A=ADATA->A[id];
313  } else {
314  DSDPSETERR2(2,"Invalid Matrix request. 0 <= %d < %d\n",id,ADATA->nnzmats);
315  }
316  DSDPFunctionReturn(0);
317 }
318 
319 #undef __FUNCT__
320 #define __FUNCT__ "DSDPBlockDataRowSparsity"
321 
330 int DSDPBlockDataRowSparsity(DSDPBlockData *ADATA,int row, int ai[], int rnnz[],int n){
331  int info,i,vari,rn;
332  DSDPFunctionBegin;
333  if (ai){
334  for (i=0; i<ADATA->nnzmats; i++){
335  vari=ADATA->nzmat[i];
336  if (ai[vari]==0){continue;}
337  info=DSDPDataMatGetRowNonzeros(ADATA->A[i],row, n, rnnz, &rn); DSDPCHKVARERR(vari,info);
338  }
339  }
340  DSDPFunctionReturn(0);
341 }
342 
343 #undef __FUNCT__
344 #define __FUNCT__ "DSDPBlockRemoveDataMatrix"
345 
352  int info,ii,k;
353  DSDPFunctionBegin;
354  for (ii=0;ii<ADATA->nnzmats;ii++){
355  if (ADATA->nzmat[ii]==vari){
356  info=DSDPDataMatDestroy(&ADATA->A[ii]);DSDPCHKVARERR(vari,info);
357  info=DSDPSetDataMatZero(&ADATA->A[ii]);DSDPCHKVARERR(vari,info);
358  for (k=ii;k<ADATA->nnzmats;k++){
359  ADATA->A[k]=ADATA->A[k+1];
360  ADATA->nzmat[k]=ADATA->nzmat[k+1];
361  }
362  ADATA->nnzmats--;
363  info=DSDPSetDataMatZero(&ADATA->A[ADATA->nnzmats]);DSDPCHKERR(info);
364  DSDPFunctionReturn(0);
365  }
366  }
367  DSDPFunctionReturn(0);
368 }
369 
370 #undef __FUNCT__
371 #define __FUNCT__ "DSDPBlockAddDataMatrix"
372 
381 int DSDPBlockAddDataMatrix(DSDPBlockData *ADATA,int vari, struct DSDPDataMat_Ops* dsdpdataops, void* data){
382  int info,ii;
383  DSDPFunctionBegin;
384  if (ADATA->nnzmats>=ADATA->maxnnzmats){
385  info=DSDPBlockDataAllocate(ADATA,2*ADATA->maxnnzmats+7);DSDPCHKERR(info);
386  }
387  ii=ADATA->nnzmats;
388  info=DSDPDataMatDestroy(&ADATA->A[ii]);DSDPCHKERR(info);
389  info=DSDPDataMatSetData(&ADATA->A[ii], dsdpdataops, data);DSDPCHKVARERR(vari,info);
390  ADATA->nzmat[ii]=vari;
391  ADATA->nnzmats++;
392  DSDPFunctionReturn(0);
393 }
394 
395 #undef __FUNCT__
396 #define __FUNCT__ "DSDPBlockSetDataMatrix"
397 
406 int DSDPBlockSetDataMatrix(DSDPBlockData *ADATA,int vari, struct DSDPDataMat_Ops* dsdpdataops, void* data){
407  int info;
408  DSDPFunctionBegin;
409  info=DSDPBlockRemoveDataMatrix(ADATA,vari);DSDPCHKERR(info);
410  info=DSDPBlockAddDataMatrix(ADATA,vari,dsdpdataops,data);DSDPCHKERR(info);
411  DSDPFunctionReturn(0);
412 }
413 
414 #undef __FUNCT__
415 #define __FUNCT__ "DSDPBlockNorm2"
416 int DSDPBlockNorm2(DSDPBlockData *ADATA, int n){
417  double fn2,tt=0;
418  int ii,info;
419  DSDPFunctionBegin;
420  for (ii=0;ii<ADATA->nnzmats;ii++){
421  info=DSDPDataMatFNorm2(ADATA->A[ii],n,&fn2); DSDPCHKERR(info);
422  tt+=fn2;
423  }
424  DSDPFunctionReturn(0);
425 }
426 
427 #undef __FUNCT__
428 #define __FUNCT__ "DSDPBlockANorm2"
429 int DSDPBlockANorm2(DSDPBlockData *ADATA, DSDPVec ANorm2, int n){
430 
431  double fn2,scl=ADATA->scl;
432  int ii,vari,info;
433 
434  DSDPFunctionBegin;
435  info=DSDPBlockNorm2(ADATA,n);DSDPCHKERR(info);
436  scl=ADATA->scl;
437  for (ii=0;ii<ADATA->nnzmats;ii++){
438  vari=ADATA->nzmat[ii];
439  info=DSDPDataMatFNorm2(ADATA->A[ii],n,&fn2); DSDPCHKVARERR(vari,info);
440  info=DSDPVecAddElement(ANorm2,vari,fn2*scl);DSDPCHKVARERR(vari,info);
441  }
442  DSDPFunctionReturn(0);
443 }
444 
445 
446 #undef __FUNCT__
447 #define __FUNCT__ "DSDPBlockView"
448 
455  int ii,kk;
456 
457  DSDPFunctionBegin;
458  for (ii=0;ii<ADATA->nnzmats;ii++){
459  kk=ADATA->nzmat[ii];
460  if (kk==0){ printf("+ C\n");}
461  else { printf(" - A[%d] y%d\n",kk,kk);}
462  }
463  printf(" = S >= 0\n");
464  DSDPFunctionReturn(0);
465 }
466 #undef __FUNCT__
467 #define __FUNCT__ "DSDPBlockView2"
468 
475  int ii,vari,info;
476 
477  DSDPFunctionBegin;
478  for (ii=0;ii<ADATA->nnzmats;ii++){
479  vari=ADATA->nzmat[ii];
480  printf("A[%d] y%d \n",vari,vari);
481  info=DSDPDataMatView(ADATA->A[ii]); DSDPCHKERR(info);
482  }
483  DSDPFunctionReturn(0);
484 }
485 
486 
487 #undef __FUNCT__
488 #define __FUNCT__ "DSDPDataMatCheck"
489 
498 int DSDPDataMatCheck(DSDPDataMat AA, SDPConeVec W, DSDPIndex IS, DSDPVMat XX){
499 
500  double *xx,ack,vAv=0,esum=0,sum,eignorm,fnorm22,dnorm,scl=1;
501  int k,n,nn,rank,info;
502 
503  DSDPFunctionBegin;
504  info=SDPConeVecGetSize(W,&n);DSDPCHKERR(info);
505 
506  info=DSDPVMatZeroEntries(XX);DSDPCHKERR(info);
507  info=DSDPDataMatGetRank(AA,&rank,n);DSDPCHKERR(info);
508  for (k=0; k<rank; k++){
509  info=DSDPDataMatGetEig(AA,k,W,IS,&ack); DSDPCHKERR(info);
510  info=SDPConeVecDot(W,W,&eignorm);DSDPCHKERR(info);
511  info=DSDPVMatAddOuterProduct(XX,scl*ack,W);DSDPCHKERR(info);
512  info=DSDPDataMatVecVec(AA,W,&sum);DSDPCHKERR(info);
513  vAv+=ack*ack*eignorm*eignorm*scl;
514  }
515  info=DSDPDataMatFNorm2(AA,n,&fnorm22); DSDPCHKERR(info);
516 
517  info=DSDPVMatScaleDiagonal(XX,0.5); DSDPCHKERR(info);
518  info=DSDPVMatGetArray(XX, &xx, &nn); DSDPCHKERR(info);
519  info=DSDPDataMatDot(AA,xx,nn,n,&esum);DSDPCHKERR(info);
520  info=DSDPVMatRestoreArray(XX, &xx, &nn); DSDPCHKERR(info);
521  info=DSDPVMatScaleDiagonal(XX,2.0); DSDPCHKERR(info);
522 
523  info=DSDPVMatGetArray(XX, &xx, &nn); DSDPCHKERR(info);
524  info=DSDPDataMatAddMultiple(AA,-1.0,xx,nn,n); DSDPCHKERR(info);
525  info=DSDPVMatRestoreArray(XX, &xx, &nn); DSDPCHKERR(info);
526  if (0==1){info=DSDPVMatView(XX);DSDPCHKERR(info);}
527  info=DSDPVMatNormF2(XX,&dnorm); DSDPCHKERR(info);
528  printf(" %4.4e, %4.4e %4.4e\n",esum,vAv,fnorm22);
529  printf(" error1: %4.4e, error2: %4.4e, error3: %4.4e\n",sqrt(dnorm),fabs(esum-vAv),fabs(fnorm22-vAv));
530  if (dnorm>1) printf("Check Add or eigs\n");
531  if (fabs(esum-vAv) > 1.0) printf("Check vAv \n");
532  if (fabs(fnorm22-vAv) > 1.0) printf("Check fnorm22\n");
533 
534  DSDPFunctionReturn(0);
535 }
536 
int DSDPDataMatDot(DSDPDataMat A, double x[], int nn, int n, double *v)
Compute inner product of data with a dense matrix.
Definition: dsdpdatamat.c:273
int DSDPBlockRemoveDataMatrix(DSDPBlockData *ADATA, int vari)
Remove a data matrix.
Definition: dsdpblock.c:351
int DSDPDataMatFNorm2(DSDPDataMat A, int n, double *fnorm2)
Compute the square of the Frobenius norm.
Definition: dsdpdatamat.c:175
struct DSDPVec_C DSDPVec
This object hold m+2 variables: a scaling of C, the y variables, and r.
Definition: dsdpvec.h:25
int DSDPBlockASum(DSDPBlockData *ADATA, double aa, DSDPVec Yk, DSDPVMat XX)
Sum the data matrices.
Definition: dsdpblock.c:20
int DSDPVMatGetArray(DSDPVMat X, double **v, int *nn)
Get the array that stores the matrix.
Definition: dsdpxmat.c:211
int DSDPBlockDataMarkNonzeroMatrices(DSDPBlockData *ADATA, int *annz)
Mark which variable in block have a data matrix.
Definition: dsdpblock.c:254
Error handling, printing, and profiling.
int DSDPBlockView(DSDPBlockData *ADATA)
Print the structure of the block.
Definition: dsdpblock.c:454
int DSDPVMatGetSize(DSDPVMat X, int *n)
Get number of rows and columns.
Definition: dsdpxmat.c:65
int DSDPVMatRestoreArray(DSDPVMat X, double **v, int *nn)
Restore the array that stores the matrix.
Definition: dsdpxmat.c:233
int DSDPDataMatGetRank(DSDPDataMat A, int *rank, int n)
Get the number of nonzero eigenvalues/eigenvectors for the matrix.
Definition: dsdpdatamat.c:129
int DSDPVMatScaleDiagonal(DSDPVMat X, double dscale)
Scaling diagonal is useful for inner products and norms.
Definition: dsdpxmat.c:147
int DSDPDataMatGetEig(DSDPDataMat A, int rr, SDPConeVec V, DSDPIndex S, double *eigenvalue)
Get an eigenvalue/vector pair.
Definition: dsdpdatamat.c:204
int DSDPDataMatCheck(DSDPDataMat AA, SDPConeVec W, DSDPIndex IS, DSDPVMat XX)
Check correctness of operations on the data.
Definition: dsdpblock.c:498
int DSDPBlockvAv(DSDPBlockData *ADATA, double aa, DSDPVec Alpha, SDPConeVec V, DSDPVec VAV)
Set VAV[i] to aa * Alpha[i] * V' A[i] V.
Definition: dsdpblock.c:84
int DSDPDataMatGetRowNonzeros(DSDPDataMat A, int nrow, int nmax, int *nz, int *nnz)
Get sparsity pattern of a row of the matrix.
Definition: dsdpdatamat.c:355
int DSDPBlockADot(DSDPBlockData *ADATA, double aa, DSDPVec Alpha, DSDPVMat X, DSDPVec AX)
Compute inner product of XX with data matrices.
Definition: dsdpblock.c:49
int DSDPBlockCountNonzeroMatrices(DSDPBlockData *ADATA, int *nzmats)
Count how many data matrices are in a block of data.
Definition: dsdpblock.c:272
int DSDPDataMatInitialize(DSDPDataMat *A)
Set pointers to NULL;.
Definition: dsdpdatamat.c:78
int DSDPBlockTakeDownData(DSDPBlockData *ADATA)
Free structures in block of data.
Definition: dsdpblock.c:182
Internal structure for data in one block of semidefintie.
Definition: dsdpsdp.h:39
Vector whose length corresponds to dimension of a block in a cone.
Definition: sdpconevec.h:13
Symmetric data matrix for one block in the semidefinite cone.
Definition: dsdpdatamat.h:15
int DSDPBlockView2(DSDPBlockData *ADATA)
Print the data.
Definition: dsdpblock.c:474
int DSDPDataMatAddMultiple(DSDPDataMat A, double ytmp, double *v, int nn, int n)
Add a multiple the data matrix to the array.
Definition: dsdpdatamat.c:402
int DSDPBlockFactorData(DSDPBlockData *ADATA, DSDPVMat X, SDPConeVec W)
Factor the data matrices.
Definition: dsdpblock.c:113
Table of function pointers that operate on the data matrix.
int DSDPDataMatDestroy(DSDPDataMat *A)
Free the data structures.
Definition: dsdpdatamat.c:444
int SDPConeVecDot(SDPConeVec V1, SDPConeVec V2, double *ans)
Inner product of two vectors.
Definition: sdpconevec.c:125
int DSDPSetDataMatZero(DSDPDataMat *A)
Make a data matrix a zero matrix.
int DSDPBlockDataAllocate(DSDPBlockData *ADATA, int nnz)
Allocate some structures.
Definition: dsdpblock.c:221
int DSDPVMatZeroEntries(DSDPVMat X)
Zero matrix.
Definition: dsdpxmat.c:125
int DSDPVMatAddOuterProduct(DSDPVMat X, double alpha, SDPConeVec V)
Add outer product of a vector to the matrix.
Definition: dsdpxmat.c:275
int DSDPBlockDataRowSparsity(DSDPBlockData *ADATA, int row, int ai[], int rnnz[], int n)
Determine sparsity pattern of data.
Definition: dsdpblock.c:330
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 DSDPDataMatVecVec(DSDPDataMat A, SDPConeVec W, double *v)
Compute w' A w.
Definition: dsdpdatamat.c:297
int DSDPVMatView(DSDPVMat X)
Print matrix.
Definition: dsdpxmat.c:107
Dense symmetric matrix for one block in the semidefinite cone.
Definition: dsdpxmat.h:17
int DSDPBlockDataDestroy(DSDPBlockData *ADATA)
Free the data matrices.
Definition: dsdpblock.c:195
int DSDPDataMatView(DSDPDataMat A)
Print matrix.
Definition: dsdpdatamat.c:423
int DSDPVMatNormF2(DSDPVMat X, double *normf2)
Compute square of Frobenius norm of matrix.
Definition: dsdpxmat.c:186
int DSDPBlockDataInitialize(DSDPBlockData *ADATA)
Set pointers to null.
Definition: dsdpblock.c:163
Internal SDPCone data structures and routines.