DSDP
sdpcompute.c
Go to the documentation of this file.
1#include "dsdpsdp.h"
2#include "dsdpsys.h"
8/*
9static int tt1=0,tt2=0;
10 DSDPEventLogBegin(tt1);
11 DSDPEventLogBegin(tt2);
12 DSDPEventLogEnd(tt2);
13 DSDPEventLogEnd(tt1);
14 if (tt1==0){DSDPEventLogRegister("SDP T1 Event",&tt1);}
15 if (tt2==0){DSDPEventLogRegister("SDP T2 Event",&tt2);}
16*/
17#undef __FUNCT__
18#define __FUNCT__ "SDPConeComputeHessian"
30int SDPConeComputeHessian( SDPCone sdpcone, double mu, DSDPSchurMat M, DSDPVec vrhs1, DSDPVec vrhs2){
31
32 int i,k,kt,kk,m,n,rank,info;
33 int ncols,ii,idA;
34 double rtemp,ack,ggamma,bmu,scl;
35 double rhs1i,rhs2i;
36 DSDPTruth method1;
37 SDPConeVec W,W2;
38 DSDPVec MRowI=sdpcone->Work, Select=sdpcone->Work2;
39 DSDPDataMat AA;
41 DSDPVMat T;
42 DSDPDataTranspose ATranspose=sdpcone->ATR;
43 SDPblk *blk=sdpcone->blk;
44 DSDPIndex IS;
45
46 /* Evaluate M */
47 DSDPFunctionBegin;
48 SDPConeValid(sdpcone);
49 info=DSDPVecGetSize(vrhs1,&m);DSDPCHKERR(info);
50
51 for (i=0; i<m; i++){ /* One row at a time */
52
53 /* Which Coluns */
54 rhs1i=0;rhs2i=0;
55 info=DSDPVecZero(MRowI);DSDPCHKERR(info);
56 info=DSDPSchurMatRowColumnScaling(M,i,Select,&ncols); DSDPCHKERR(info);
57 if (ncols==0){continue;}
58
59 for (kt=0; kt<ATranspose.nnzblocks[i]; kt++){ /* Loop over all blocks */
60 kk=ATranspose.nzblocks[i][kt];
61 idA=ATranspose.idA[i][kt];
62 info=DSDPBlockGetMatrix(&blk[kk].ADATA,idA,&ii,&scl,&AA);DSDPCHKBLOCKERR(kk,info);
63 if (ii!=i){DSDPSETERR2(8,"Data Transpose Error: var %d does not equal %d.\n",i,ii);}
64 info = DSDPDataMatGetRank(AA,&rank,blk[kk].n);DSDPCHKBLOCKERR(kk,info);
65 if (rank==0) continue;
66
67 T=blk[kk].T; S=blk[kk].S; W=blk[kk].W; W2=blk[kk].W2;
68 n=blk[kk].n; ggamma=blk[kk].gammamu; bmu=blk[kk].bmu; IS=blk[kk].IS;
69
70 method1=DSDP_TRUE; /* Simple heuristic */
71 if (rank==1) method1=DSDP_FALSE;
72 if (rank==2 && ncols<=n) method1=DSDP_FALSE;
73 if (rank*rank*ncols<=n+1)method1=DSDP_FALSE;
74 if (ncols*blk[kk].nnz < n*n/10) method1=DSDP_FALSE;
75 if (ncols==1 && i==m-1)method1=DSDP_FALSE;
76 if (n<5) method1=DSDP_TRUE;
77 if (0==1) method1=DSDP_FALSE;
78 if (method1==DSDP_TRUE){info=DSDPVMatZeroEntries(T);DSDPCHKBLOCKERR(kk,info);}
79 for (k=0; k<rank; k++){
80
81 info=DSDPDataMatGetEig(AA,k,W,IS,&ack); DSDPCHKBLOCKERR(kk,info);
82 if (ack==0.0) continue;
83 ack*=scl;
84 info=DSDPDualMatInverseMultiply(S,IS,W,W2);DSDPCHKBLOCKERR(kk,info);
85
86 /* RHS terms */
87 info = SDPConeVecDot(W,W2,&rtemp); DSDPCHKBLOCKERR(kk,info);
88 if (rtemp==0.0) continue;
89 rhs1i+=rtemp*ack*bmu; rhs2i+=rtemp*ack*ggamma*mu;
90 ack*=(ggamma+bmu);
91
92 if (method1==DSDP_TRUE){
93 info=DSDPVMatAddOuterProduct(T,ack*mu,W2);DSDPCHKBLOCKERR(kk,info);
94 } else {
95 info=DSDPBlockvAv(&blk[kk].ADATA,ack*mu,Select,W2,MRowI);DSDPCHKBLOCKERR(kk,info);
96 } /* End row computations for rank kk of block kk */
97
98 } /* End row computations for all of block kk */
99
100 if (method1==DSDP_TRUE){
101 info=DSDPBlockADot(&blk[kk].ADATA,1.0,Select,T,MRowI);DSDPCHKBLOCKERR(kk,info);
102 } /* End row computations for all of block ll */
103 } /* End row computations for all blocks */
104 info=DSDPVecAddElement(vrhs1,i,rhs1i);DSDPCHKERR(info);
105 info=DSDPVecAddElement(vrhs2,i,rhs2i);DSDPCHKERR(info);
106 info=DSDPSchurMatAddRow(M,i,1.0,MRowI);DSDPCHKERR(info);
107 }
108
109 DSDPFunctionReturn(0);
110}
111
112#undef __FUNCT__
113#define __FUNCT__ "SDPConeComputeRHS"
125int SDPConeComputeRHS( SDPCone sdpcone, int blockj, double mu, DSDPVec vrow, DSDPVec vrhs1, DSDPVec vrhs2){
126
127 int info,i,ii,k,rank,nnzmats;
128 double dtmp,dyiscale=1,ack,scl,rtemp;
129 SDPblk *sdp=&sdpcone->blk[blockj];
130 SDPConeVec W=sdp->W,W2=sdp->W2;
131 DSDPDataMat AA;
132 DSDPVMat T=sdp->T;
133 DSDPDualMat S=sdp->S;
134 DSDPTruth method1;
135 DSDPIndex IS=sdp->IS;
136
137 DSDPFunctionBegin;
138
139 info=SDPConeCheckJ(sdpcone,blockj);DSDPCHKERR(info);
140 method1=DSDP_TRUE;
141 method1=DSDP_FALSE;
142
143 if (method1==DSDP_TRUE){
144 info=DSDPBlockCountNonzeroMatrices(&sdp->ADATA,&nnzmats);DSDPCHKERR(info);
145 for (i=0;i<nnzmats;i++){
146 info=DSDPBlockGetMatrix(&sdp->ADATA,i,&ii,&scl,&AA);DSDPCHKERR(info);
147 info=DSDPVecGetElement(vrow,ii,&dyiscale);DSDPCHKVARERR(ii,info);
148 if (dyiscale==0) continue;
149 info=DSDPDataMatGetRank(AA,&rank,sdp->n);DSDPCHKVARERR(ii,info);
150 for (k=0; k<rank; k++){
151 info=DSDPDataMatGetEig(AA,k,W,IS,&ack); DSDPCHKVARERR(ii,info);
152 if (ack==0) continue;
153 info=DSDPDualMatInverseMultiply(S,IS,W,W2);DSDPCHKVARERR(ii,info);
154 info=SDPConeVecDot(W,W2,&rtemp); DSDPCHKVARERR(ii,info);
155 dtmp=rtemp*ack*mu*dyiscale*scl;
156 info=DSDPVecAddElement(vrhs2,ii,dtmp);DSDPCHKVARERR(ii,info);
157 }
158 }
159
160 } else {
161 info=DSDPVMatZeroEntries(T);DSDPCHKERR(info);
162 info=DSDPDualMatInverseAdd(S,mu,T);DSDPCHKERR(info);
163 info=DSDPBlockADot(&sdp->ADATA,1.0,vrow,T,vrhs2);DSDPCHKERR(info);
164 }
165
166 DSDPFunctionReturn(0);
167}
168
169#undef __FUNCT__
170#define __FUNCT__ "SDPConeMultiply"
182int SDPConeMultiply(SDPCone sdpcone, int blockj, double mu, DSDPVec vrow, DSDPVec vin, DSDPVec vout){
183
184 int info,i,ii,k,rank,nnzmats;
185 double dd2,dtmp,dyiscale,ack,scl,vv;
186 SDPblk *sdp=&sdpcone->blk[blockj];
187 SDPConeVec W=sdp->W,W2=sdp->W2;
188 DSDPDataMat AA;
189 DSDPVMat T=sdp->T;
190 DSDPDSMat DS=sdp->DS;
191 DSDPIndex IS=sdp->IS;
192 DSDPDualMat S=sdp->S;
193
194 DSDPFunctionBegin;
195 info=SDPConeCheckJ(sdpcone,blockj);DSDPCHKERR(info);
196 info=DSDPVMatZeroEntries(T); DSDPCHKERR(info);
197 info=DSDPBlockASum(&sdp->ADATA,-1.0,vin,T); DSDPCHKERR(info);
198 info=DSDPDSMatSetArray(DS,T); DSDPCHKERR(info);
199 info=DSDPBlockCountNonzeroMatrices(&sdp->ADATA,&nnzmats);DSDPCHKERR(info);
200 for (i=0;i<nnzmats;i++){
201 info=DSDPBlockGetMatrix(&sdp->ADATA,i,&ii,&scl,&AA);DSDPCHKERR(info);
202 info=DSDPVecGetElement(vrow,ii,&dyiscale);DSDPCHKVARERR(ii,info);
203 if (dyiscale==0) continue;
204
205 info=DSDPDataMatGetRank(AA,&rank,sdp->n);DSDPCHKVARERR(ii,info);
206 for (dd2=0,k=0; k<rank; k++){
207 info=DSDPDataMatGetEig(AA,k,W,IS,&ack); DSDPCHKVARERR(ii,info);
208 if (ack==0) continue;
209 info=DSDPDualMatInverseMultiply(S,IS,W,W2);DSDPCHKVARERR(ii,info);
210 info=DSDPDSMatVecVec(DS,W2,&vv);DSDPCHKVARERR(ii,info);
211 dd2+=vv*ack;
212 }
213 dtmp = dd2 * dyiscale *mu *scl;
214 info=DSDPVecAddElement(vout,ii,dtmp);DSDPCHKVARERR(ii,info);
215 }
216
217 DSDPFunctionReturn(0);
218}
219
220
221
222#undef __FUNCT__
223#define __FUNCT__ "SDPConeComputeXX"
235int SDPConeComputeXX(SDPCone sdpcone, int blockj, DSDPVec DY, double mu,DSDPDualMat S, DSDPVMat X){
236
237 int info, i, ii,k, rank, n, nnzmats;
238 double dtmp,dyiscale,ack,scl;
239 SDPblk *sdp=&sdpcone->blk[blockj];
240 SDPConeVec W=sdp->W,W2=sdp->W2;
241 DSDPDataMat AA;
242 DSDPIndex IS=sdp->IS;
243
244 DSDPFunctionBegin;
245 info=SDPConeCheckJ(sdpcone,blockj);DSDPCHKERR(info);
246 mu=mu*sdp->gammamu; n=sdp->n;
247 info=DSDPVMatZeroEntries(X);DSDPCHKERR(info);
248 info=DSDPBlockCountNonzeroMatrices(&sdp->ADATA,&nnzmats);DSDPCHKERR(info);
249 for (i=0;i<nnzmats;i++){
250 info=DSDPBlockGetMatrix(&sdp->ADATA,i,&ii,&scl,&AA);DSDPCHKVARERR(ii,info);
251 info=DSDPVecGetElement(DY,ii,&dyiscale);DSDPCHKVARERR(ii,info);
252 if (dyiscale==0) continue;
253 info=DSDPDataMatGetRank(AA,&rank,n);DSDPCHKVARERR(ii,info);
254 for (k=0; k<rank; k++){
255 info = DSDPDataMatGetEig(AA,k,W,IS,&ack);DSDPCHKVARERR(ii,info);
256 if (ack==0) continue;
257 info=DSDPDualMatInverseMultiply(S,IS,W,W2);DSDPCHKVARERR(ii,info);
258 dtmp = ack * dyiscale * mu * scl;
259 info=DSDPVMatAddOuterProduct(X,dtmp,W2);DSDPCHKVARERR(ii,info);
260 }
261 }
262
263 info=DSDPDualMatInverseAdd(S,mu,X);DSDPCHKERR(info);
264 DSDPFunctionReturn(0);
265}
266
int SDPConeCheckJ(SDPCone sdpcone, int blockj)
Check validity of parameter.
Definition dsdpadddata.c:31
DSDPTruth
Boolean variables.
@ DSDP_FALSE
@ DSDP_TRUE
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 DSDPBlockASum(DSDPBlockData *ADATA, double aa, DSDPVec Yk, DSDPVMat XX)
Sum the data matrices.
Definition dsdpblock.c:20
int DSDPBlockCountNonzeroMatrices(DSDPBlockData *ADATA, int *nzmats)
Count how many data matrices are in a block of data.
Definition dsdpblock.c:272
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 DSDPBlockADot(DSDPBlockData *ADATA, double aa, DSDPVec Alpha, DSDPVMat X, DSDPVec AX)
Compute inner product of XX with data matrices.
Definition dsdpblock.c:49
int DSDPDataMatGetEig(DSDPDataMat A, int rr, SDPConeVec V, DSDPIndex S, double *eigenvalue)
Get an eigenvalue/vector pair.
int DSDPDataMatGetRank(DSDPDataMat A, int *rank, int n)
Get the number of nonzero eigenvalues/eigenvectors for the matrix.
int DSDPDSMatVecVec(DSDPDSMat A, SDPConeVec X, double *vAv)
Compute the product x' A x.
Definition dsdpdsmat.c:181
int DSDPDSMatSetArray(DSDPDSMat A, DSDPVMat T)
Set values into the matrix.
Definition dsdpdsmat.c:130
int DSDPDualMatInverseAdd(DSDPDualMat S, double alpha, DSDPVMat T)
Add a multiple of the inverse to T.
int DSDPDualMatInverseMultiply(DSDPDualMat S, DSDPIndex IS, SDPConeVec B, SDPConeVec X)
Multiply the inverse by a vector or solve the system of equations.
int DSDPSchurMatRowColumnScaling(DSDPSchurMat, int, DSDPVec, int *)
Get the scaling and nonzero pattern of each column in this row of the matrix.
int DSDPSchurMatAddRow(DSDPSchurMat, int, double, DSDPVec)
Add elements to a row of the Schur matrix.
Internal SDPCone data structures and routines.
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 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 SDPConeComputeHessian(SDPCone sdpcone, double mu, DSDPSchurMat M, DSDPVec vrhs1, DSDPVec vrhs2)
Compute the Hessian to the barrier term.
Definition sdpcompute.c:30
int SDPConeComputeXX(SDPCone sdpcone, int blockj, DSDPVec DY, double mu, DSDPDualMat S, DSDPVMat X)
Compute X.
Definition sdpcompute.c:235
int SDPConeMultiply(SDPCone sdpcone, int blockj, double mu, DSDPVec vrow, DSDPVec vin, DSDPVec vout)
Compute the gradient to the barrier term.
Definition sdpcompute.c:182
int SDPConeComputeRHS(SDPCone sdpcone, int blockj, double mu, DSDPVec vrow, DSDPVec vrhs1, DSDPVec vrhs2)
Compute the gradient to the barrier term.
Definition sdpcompute.c:125
int SDPConeVecDot(SDPConeVec V1, SDPConeVec V2, double *ans)
Inner product of two vectors.
Definition sdpconevec.c:125
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.
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
Internal structure for semidefinite cone.
Definition dsdpsdp.h:80
Internal structure for block of semidefinite cone.
Definition dsdpsdp.h:52