sdpcompute.c
Go to the documentation of this file.
#include "dsdpsdp.h"
#include "dsdpsys.h"
/*
static int tt1=0,tt2=0;
DSDPEventLogBegin(tt1);
DSDPEventLogBegin(tt2);
DSDPEventLogEnd(tt2);
DSDPEventLogEnd(tt1);
if (tt1==0){DSDPEventLogRegister("SDP T1 Event",&tt1);}
if (tt2==0){DSDPEventLogRegister("SDP T2 Event",&tt2);}
*/
#undef __FUNCT__
#define __FUNCT__ "SDPConeComputeHessian"
int SDPConeComputeHessian( SDPCone sdpcone, double mu, DSDPSchurMat M, DSDPVec vrhs1, DSDPVec vrhs2){
int i,k,kt,kk,m,n,rank,info;
int ncols,ii,idA;
double rtemp,ack,ggamma,bmu,scl;
double rhs1i,rhs2i;
DSDPTruth method1;
SDPConeVec W,W2;
DSDPVec MRowI=sdpcone->Work, Select=sdpcone->Work2;
DSDPDataMat AA;
DSDPDualMat S;
DSDPVMat T;
DSDPDataTranspose ATranspose=sdpcone->ATR;
SDPblk *blk=sdpcone->blk;
DSDPIndex IS;
/* Evaluate M */
DSDPFunctionBegin;
SDPConeValid(sdpcone);
info=DSDPVecGetSize(vrhs1,&m);DSDPCHKERR(info);
for (i=0; i<m; i++){ /* One row at a time */
/* Which Coluns */
rhs1i=0;rhs2i=0;
info=DSDPVecZero(MRowI);DSDPCHKERR(info);
info=DSDPSchurMatRowColumnScaling(M,i,Select,&ncols); DSDPCHKERR(info);
if (ncols==0){continue;}
for (kt=0; kt<ATranspose.nnzblocks[i]; kt++){ /* Loop over all blocks */
kk=ATranspose.nzblocks[i][kt];
idA=ATranspose.idA[i][kt];
info=DSDPBlockGetMatrix(&blk[kk].ADATA,idA,&ii,&scl,&AA);DSDPCHKBLOCKERR(kk,info);
if (ii!=i){DSDPSETERR2(8,"Data Transpose Error: var %d does not equal %d.\n",i,ii);}
info = DSDPDataMatGetRank(AA,&rank,blk[kk].n);DSDPCHKBLOCKERR(kk,info);
if (rank==0) continue;
T=blk[kk].T; S=blk[kk].S; W=blk[kk].W; W2=blk[kk].W2;
n=blk[kk].n; ggamma=blk[kk].gammamu; bmu=blk[kk].bmu; IS=blk[kk].IS;
method1=DSDP_TRUE; /* Simple heuristic */
if (rank==1) method1=DSDP_FALSE;
if (rank==2 && ncols<=n) method1=DSDP_FALSE;
if (rank*rank*ncols<=n+1)method1=DSDP_FALSE;
if (ncols*blk[kk].nnz < n*n/10) method1=DSDP_FALSE;
if (ncols==1 && i==m-1)method1=DSDP_FALSE;
if (n<5) method1=DSDP_TRUE;
if (0==1) method1=DSDP_FALSE;
if (method1==DSDP_TRUE){info=DSDPVMatZeroEntries(T);DSDPCHKBLOCKERR(kk,info);}
for (k=0; k<rank; k++){
info=DSDPDataMatGetEig(AA,k,W,IS,&ack); DSDPCHKBLOCKERR(kk,info);
if (ack==0.0) continue;
ack*=scl;
info=DSDPDualMatInverseMultiply(S,IS,W,W2);DSDPCHKBLOCKERR(kk,info);
/* RHS terms */
info = SDPConeVecDot(W,W2,&rtemp); DSDPCHKBLOCKERR(kk,info);
if (rtemp==0.0) continue;
rhs1i+=rtemp*ack*bmu; rhs2i+=rtemp*ack*ggamma*mu;
ack*=(ggamma+bmu);
if (method1==DSDP_TRUE){
info=DSDPVMatAddOuterProduct(T,ack*mu,W2);DSDPCHKBLOCKERR(kk,info);
} else {
info=DSDPBlockvAv(&blk[kk].ADATA,ack*mu,Select,W2,MRowI);DSDPCHKBLOCKERR(kk,info);
} /* End row computations for rank kk of block kk */
} /* End row computations for all of block kk */
if (method1==DSDP_TRUE){
info=DSDPBlockADot(&blk[kk].ADATA,1.0,Select,T,MRowI);DSDPCHKBLOCKERR(kk,info);
} /* End row computations for all of block ll */
} /* End row computations for all blocks */
info=DSDPVecAddElement(vrhs1,i,rhs1i);DSDPCHKERR(info);
info=DSDPVecAddElement(vrhs2,i,rhs2i);DSDPCHKERR(info);
info=DSDPSchurMatAddRow(M,i,1.0,MRowI);DSDPCHKERR(info);
}
DSDPFunctionReturn(0);
}
#undef __FUNCT__
#define __FUNCT__ "SDPConeComputeRHS"
int SDPConeComputeRHS( SDPCone sdpcone, int blockj, double mu, DSDPVec vrow, DSDPVec vrhs1, DSDPVec vrhs2){
int info,i,ii,k,rank,nnzmats;
double dtmp,dyiscale=1,ack,scl,rtemp;
SDPblk *sdp=&sdpcone->blk[blockj];
SDPConeVec W=sdp->W,W2=sdp->W2;
DSDPDataMat AA;
DSDPVMat T=sdp->T;
DSDPDualMat S=sdp->S;
DSDPTruth method1;
DSDPIndex IS=sdp->IS;
DSDPFunctionBegin;
info=SDPConeCheckJ(sdpcone,blockj);DSDPCHKERR(info);
method1=DSDP_TRUE;
method1=DSDP_FALSE;
if (method1==DSDP_TRUE){
info=DSDPBlockCountNonzeroMatrices(&sdp->ADATA,&nnzmats);DSDPCHKERR(info);
for (i=0;i<nnzmats;i++){
info=DSDPBlockGetMatrix(&sdp->ADATA,i,&ii,&scl,&AA);DSDPCHKERR(info);
info=DSDPVecGetElement(vrow,ii,&dyiscale);DSDPCHKVARERR(ii,info);
if (dyiscale==0) continue;
info=DSDPDataMatGetRank(AA,&rank,sdp->n);DSDPCHKVARERR(ii,info);
for (k=0; k<rank; k++){
info=DSDPDataMatGetEig(AA,k,W,IS,&ack); DSDPCHKVARERR(ii,info);
if (ack==0) continue;
info=DSDPDualMatInverseMultiply(S,IS,W,W2);DSDPCHKVARERR(ii,info);
info=SDPConeVecDot(W,W2,&rtemp); DSDPCHKVARERR(ii,info);
dtmp=rtemp*ack*mu*dyiscale*scl;
info=DSDPVecAddElement(vrhs2,ii,dtmp);DSDPCHKVARERR(ii,info);
}
}
} else {
info=DSDPVMatZeroEntries(T);DSDPCHKERR(info);
info=DSDPDualMatInverseAdd(S,mu,T);DSDPCHKERR(info);
info=DSDPBlockADot(&sdp->ADATA,1.0,vrow,T,vrhs2);DSDPCHKERR(info);
}
DSDPFunctionReturn(0);
}
#undef __FUNCT__
#define __FUNCT__ "SDPConeMultiply"
int SDPConeMultiply(SDPCone sdpcone, int blockj, double mu, DSDPVec vrow, DSDPVec vin, DSDPVec vout){
int info,i,ii,k,rank,nnzmats;
double dd2,dtmp,dyiscale,ack,scl,vv;
SDPblk *sdp=&sdpcone->blk[blockj];
SDPConeVec W=sdp->W,W2=sdp->W2;
DSDPDataMat AA;
DSDPVMat T=sdp->T;
DSDPDSMat DS=sdp->DS;
DSDPIndex IS=sdp->IS;
DSDPDualMat S=sdp->S;
DSDPFunctionBegin;
info=SDPConeCheckJ(sdpcone,blockj);DSDPCHKERR(info);
info=DSDPVMatZeroEntries(T); DSDPCHKERR(info);
info=DSDPBlockASum(&sdp->ADATA,-1.0,vin,T); DSDPCHKERR(info);
info=DSDPDSMatSetArray(DS,T); DSDPCHKERR(info);
info=DSDPBlockCountNonzeroMatrices(&sdp->ADATA,&nnzmats);DSDPCHKERR(info);
for (i=0;i<nnzmats;i++){
info=DSDPBlockGetMatrix(&sdp->ADATA,i,&ii,&scl,&AA);DSDPCHKERR(info);
info=DSDPVecGetElement(vrow,ii,&dyiscale);DSDPCHKVARERR(ii,info);
if (dyiscale==0) continue;
info=DSDPDataMatGetRank(AA,&rank,sdp->n);DSDPCHKVARERR(ii,info);
for (dd2=0,k=0; k<rank; k++){
info=DSDPDataMatGetEig(AA,k,W,IS,&ack); DSDPCHKVARERR(ii,info);
if (ack==0) continue;
info=DSDPDualMatInverseMultiply(S,IS,W,W2);DSDPCHKVARERR(ii,info);
info=DSDPDSMatVecVec(DS,W2,&vv);DSDPCHKVARERR(ii,info);
dd2+=vv*ack;
}
dtmp = dd2 * dyiscale *mu *scl;
info=DSDPVecAddElement(vout,ii,dtmp);DSDPCHKVARERR(ii,info);
}
DSDPFunctionReturn(0);
}
#undef __FUNCT__
#define __FUNCT__ "SDPConeComputeXX"
int SDPConeComputeXX(SDPCone sdpcone, int blockj, DSDPVec DY, double mu,DSDPDualMat S, DSDPVMat X){
int info, i, ii,k, rank, n, nnzmats;
double dtmp,dyiscale,ack,scl;
SDPblk *sdp=&sdpcone->blk[blockj];
SDPConeVec W=sdp->W,W2=sdp->W2;
DSDPDataMat AA;
DSDPIndex IS=sdp->IS;
DSDPFunctionBegin;
info=SDPConeCheckJ(sdpcone,blockj);DSDPCHKERR(info);
mu=mu*sdp->gammamu; n=sdp->n;
info=DSDPVMatZeroEntries(X);DSDPCHKERR(info);
info=DSDPBlockCountNonzeroMatrices(&sdp->ADATA,&nnzmats);DSDPCHKERR(info);
for (i=0;i<nnzmats;i++){
info=DSDPBlockGetMatrix(&sdp->ADATA,i,&ii,&scl,&AA);DSDPCHKVARERR(ii,info);
info=DSDPVecGetElement(DY,ii,&dyiscale);DSDPCHKVARERR(ii,info);
if (dyiscale==0) continue;
info=DSDPDataMatGetRank(AA,&rank,n);DSDPCHKVARERR(ii,info);
for (k=0; k<rank; k++){
info = DSDPDataMatGetEig(AA,k,W,IS,&ack);DSDPCHKVARERR(ii,info);
if (ack==0) continue;
info=DSDPDualMatInverseMultiply(S,IS,W,W2);DSDPCHKVARERR(ii,info);
dtmp = ack * dyiscale * mu * scl;
info=DSDPVMatAddOuterProduct(X,dtmp,W2);DSDPCHKVARERR(ii,info);
}
}
info=DSDPDualMatInverseAdd(S,mu,X);DSDPCHKERR(info);
DSDPFunctionReturn(0);
}
SDPConeCheckJ
int SDPConeCheckJ(SDPCone sdpcone, int blockj)
Check validity of parameter.
Definition dsdpadddata.c:31
DSDPTruth
DSDPTruth
Boolean variables.
Definition dsdpbasictypes.h:19
DSDP_FALSE
@ DSDP_FALSE
Definition dsdpbasictypes.h:19
DSDP_TRUE
@ DSDP_TRUE
Definition dsdpbasictypes.h:19
DSDPBlockGetMatrix
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
DSDPBlockASum
int DSDPBlockASum(DSDPBlockData *ADATA, double aa, DSDPVec Yk, DSDPVMat XX)
Sum the data matrices.
Definition dsdpblock.c:20
DSDPBlockCountNonzeroMatrices
int DSDPBlockCountNonzeroMatrices(DSDPBlockData *ADATA, int *nzmats)
Count how many data matrices are in a block of data.
Definition dsdpblock.c:272
DSDPBlockvAv
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
DSDPBlockADot
int DSDPBlockADot(DSDPBlockData *ADATA, double aa, DSDPVec Alpha, DSDPVMat X, DSDPVec AX)
Compute inner product of XX with data matrices.
Definition dsdpblock.c:49
DSDPDataMatGetEig
int DSDPDataMatGetEig(DSDPDataMat A, int rr, SDPConeVec V, DSDPIndex S, double *eigenvalue)
Get an eigenvalue/vector pair.
Definition dsdpdatamat.c:204
DSDPDataMatGetRank
int DSDPDataMatGetRank(DSDPDataMat A, int *rank, int n)
Get the number of nonzero eigenvalues/eigenvectors for the matrix.
Definition dsdpdatamat.c:129
DSDPDSMatVecVec
int DSDPDSMatVecVec(DSDPDSMat A, SDPConeVec X, double *vAv)
Compute the product x' A x.
Definition dsdpdsmat.c:181
DSDPDSMatSetArray
int DSDPDSMatSetArray(DSDPDSMat A, DSDPVMat T)
Set values into the matrix.
Definition dsdpdsmat.c:130
DSDPDualMatInverseAdd
int DSDPDualMatInverseAdd(DSDPDualMat S, double alpha, DSDPVMat T)
Add a multiple of the inverse to T.
Definition dsdpdualmat.c:209
DSDPDualMatInverseMultiply
int DSDPDualMatInverseMultiply(DSDPDualMat S, DSDPIndex IS, SDPConeVec B, SDPConeVec X)
Multiply the inverse by a vector or solve the system of equations.
Definition dsdpdualmat.c:236
DSDPSchurMatRowColumnScaling
int DSDPSchurMatRowColumnScaling(DSDPSchurMat, int, DSDPVec, int *)
Get the scaling and nonzero pattern of each column in this row of the matrix.
Definition dsdpschurmatadd.c:33
DSDPSchurMatAddRow
int DSDPSchurMatAddRow(DSDPSchurMat, int, double, DSDPVec)
Add elements to a row of the Schur matrix.
Definition dsdpschurmatadd.c:76
dsdpsdp.h
Internal SDPCone data structures and routines.
dsdpsys.h
Error handling, printing, and profiling.
DSDPVec
struct DSDPVec_C DSDPVec
This object hold m+2 variables: a scaling of C, the y variables, and r.
Definition dsdpvec.h:25
DSDPVMatZeroEntries
int DSDPVMatZeroEntries(DSDPVMat X)
Zero matrix.
Definition dsdpxmat.c:125
DSDPVMatAddOuterProduct
int DSDPVMatAddOuterProduct(DSDPVMat X, double alpha, SDPConeVec V)
Add outer product of a vector to the matrix.
Definition dsdpxmat.c:275
SDPConeComputeHessian
int SDPConeComputeHessian(SDPCone sdpcone, double mu, DSDPSchurMat M, DSDPVec vrhs1, DSDPVec vrhs2)
Compute the Hessian to the barrier term.
Definition sdpcompute.c:30
SDPConeComputeXX
int SDPConeComputeXX(SDPCone sdpcone, int blockj, DSDPVec DY, double mu, DSDPDualMat S, DSDPVMat X)
Compute X.
Definition sdpcompute.c:235
SDPConeMultiply
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
SDPConeComputeRHS
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
SDPConeVecDot
int SDPConeVecDot(SDPConeVec V1, SDPConeVec V2, double *ans)
Inner product of two vectors.
Definition sdpconevec.c:125
DSDPDSMat_C
Symmetric Delta S matrix for one block in the semidefinite cone.
Definition dsdpdsmat.h:23
DSDPDataMat_C
Symmetric data matrix for one block in the semidefinite cone.
Definition dsdpdatamat.h:15
DSDPDataTranspose
Internal structure for transpose of data.
Definition dsdpsdp.h:25
DSDPDualMat_C
Represents an S matrix for one block in the semidefinite cone.
Definition dsdpdualmat.h:18
DSDPSchurMat_C
Schur complement matrix whose solution is the Newton direction.
Definition dsdpschurmat.h:35
DSDPVMat_C
Dense symmetric matrix for one block in the semidefinite cone.
Definition dsdpxmat.h:17
SDPConeVec_C
Vector whose length corresponds to dimension of a block in a cone.
Definition sdpconevec.h:13
SDPCone_C
Internal structure for semidefinite cone.
Definition dsdpsdp.h:80
SDPblk
Internal structure for block of semidefinite cone.
Definition dsdpsdp.h:52