sdpcone.c
Go to the documentation of this file.
#include "dsdpsdp.h"
#include "dsdpsys.h"
#undef __FUNCT__
#define __FUNCT__ "SDPConeComputeSS"
int SDPConeComputeSS(SDPCone sdpcone, int blockj, DSDPVec Y, DSDPVMat SS){
int info;
DSDPFunctionBegin;
info=DSDPVMatZeroEntries(SS); DSDPCHKBLOCKERR(blockj,info);
info=DSDPBlockASum(&sdpcone->blk[blockj].ADATA,1,Y,SS); DSDPCHKBLOCKERR(blockj,info);
DSDPFunctionReturn(0);
}
#undef __FUNCT__
#define __FUNCT__ "SDPConeComputeS"
int SDPConeComputeS(SDPCone sdpcone, int blockj, double cc,double y[], int nvars, double r, int n, double s[], int nn){
int i,info;
char UPLQ;
DSDPVMat T;
DSDPVec Y=sdpcone->Work;
DSDPFunctionBegin;
info=SDPConeCheckN(sdpcone,blockj,n);DSDPCHKBLOCKERR(blockj,info);
info=SDPConeCheckM(sdpcone,nvars);DSDPCHKERR(info);
if (n<1){DSDPFunctionReturn(0);}
info=DSDPVecSetC(Y,-1.0*cc);
info=DSDPVecSetR(Y,-r);
for (i=0;i<nvars;i++){info=DSDPVecSetElement(Y,i+1,y[i]);}
info=SDPConeGetStorageFormat(sdpcone,blockj,&UPLQ);DSDPCHKBLOCKERR(blockj,info);
info=DSDPMakeVMatWithArray(UPLQ,s,nn,n,&T);DSDPCHKBLOCKERR(blockj,info);
info=SDPConeComputeSS(sdpcone,blockj,Y,T);DSDPCHKBLOCKERR(blockj,info);
info=DSDPVMatDestroy(&T);DSDPCHKBLOCKERR(blockj,info);
DSDPFunctionReturn(0);
}
#undef __FUNCT__
#define __FUNCT__ "SDPConeAddADotX"
int SDPConeAddADotX(SDPCone sdpcone, int blockj, double alpha,double x[], int nn, double adotx[], int m){
int info,n;
char UPLQ;
SDPblk *blk=sdpcone->blk;
double scl=blk[blockj].ADATA.scl;
DSDPVec ADOTX,YW2;
DSDPVMat T;
DSDPFunctionBegin;
info=SDPConeCheckJ(sdpcone,blockj);DSDPCHKERR(info);
info=SDPConeCheckM(sdpcone,m-2);DSDPCHKERR(info);
YW2=sdpcone->Work2;
info=DSDPVecSet(alpha,YW2);DSDPCHKBLOCKERR(blockj,info);
info=SDPConeGetBlockSize(sdpcone,blockj,&n);DSDPCHKBLOCKERR(blockj,info);
if (n<1){DSDPFunctionReturn(0);}
info=DSDPVecCreateWArray(&ADOTX,adotx,m);
info=SDPConeGetStorageFormat(sdpcone,blockj,&UPLQ);DSDPCHKBLOCKERR(blockj,info);
info=DSDPMakeVMatWithArray(UPLQ,x,nn,n,&T);DSDPCHKBLOCKERR(blockj,info);
info=DSDPBlockADot(&blk[blockj].ADATA,1.0/scl,YW2,T,ADOTX);DSDPCHKBLOCKERR(blockj,info);
info=DSDPVMatDestroy(&T);DSDPCHKBLOCKERR(blockj,info);
DSDPFunctionReturn(0);
}
#undef __FUNCT__
#define __FUNCT__ "SDPConeComputeXDot"
int SDPConeComputeXDot(SDPCone sdpcone, int blockj, DSDPVec Y,DSDPVMat X, DSDPVec AX, double* xtrace,double *xnorm, double *tracexs){
int info;
SDPblk *blk=sdpcone->blk;
DSDPVec YW2=sdpcone->Work2;
double one=1.0,scl=blk[blockj].ADATA.scl;
DSDPFunctionBegin;
info=DSDPVecZero(YW2);DSDPCHKBLOCKERR(blockj,info);
info=DSDPBlockADot(&blk[blockj].ADATA,-1.0/scl,Y,X,YW2);DSDPCHKBLOCKERR(blockj,info);
info=DSDPVecGetR(YW2,xtrace);DSDPCHKBLOCKERR(blockj,info);
info=DSDPVecSum(YW2,tracexs);DSDPCHKBLOCKERR(blockj,info);
info=DSDPVMatNormF2(X,xnorm);DSDPCHKBLOCKERR(blockj,info);
info=DSDPVecSet(one,YW2);DSDPCHKBLOCKERR(blockj,info);
info=DSDPBlockADot(&blk[blockj].ADATA,1.0/scl,YW2,X,AX);DSDPCHKBLOCKERR(blockj,info);
DSDPFunctionReturn(0);
}
#undef __FUNCT__
#define __FUNCT__ "SDPConeComputeX3"
int SDPConeComputeX3(SDPCone sdpcone, int blockj, double mu, DSDPVec Y,DSDPVec DY,DSDPVMat X){
int info;
double xshift=1e-12,xscale=1e-12;
SDPblk *blk=sdpcone->blk;
DSDPTruth psdefinite1=DSDP_FALSE,psdefinite2=DSDP_FALSE,full;
DSDPDualMat SS;
DSDPFunctionBegin;
SS=blk[blockj].SS;
info=SDPConeComputeSS(sdpcone,blockj,Y,X);DSDPCHKBLOCKERR(blockj,info);
info=DSDPDualMatSetArray(SS,X); DSDPCHKBLOCKERR(blockj,info);
info=DSDPDualMatCholeskyFactor(SS,&psdefinite1); DSDPCHKBLOCKERR(blockj,info);
if (psdefinite1 == DSDP_FALSE){
DSDPLogInfo(0,2,"Primal SDP Block %2.0f not PSD\n",blockj);
}
info=DSDPDualMatInvert(SS);DSDPCHKBLOCKERR(blockj,info);
info=SDPConeComputeXX(sdpcone,blockj,DY,mu,SS,X);DSDPCHKBLOCKERR(blockj,info);
info=DSDPDualMatIsFull(SS,&full);DSDPCHKBLOCKERR(blockj,info);
psdefinite2=DSDP_FALSE;
while (full==DSDP_TRUE && psdefinite2==DSDP_FALSE && xscale<2e-1){
info=DSDPVMatShiftDiagonal(X,xshift); DSDPCHKBLOCKERR(blockj,info);
info=DSDPVMatScaleDiagonal(X,1.0+xscale); DSDPCHKBLOCKERR(blockj,info);
DSDPLogInfo(0,10,"VMat: shift diagonal: %4.2e, scale diagonal: %4.2e.\n",xshift,1+xscale);
info=DSDPDualMatSetArray(SS,X); DSDPCHKBLOCKERR(blockj,info);
info=DSDPDualMatCholeskyFactor(SS,&psdefinite2); DSDPCHKBLOCKERR(blockj,info);
xshift*=10;xscale*=10;
}
if (full==DSDP_FALSE){
xshift=1e-12,xscale=1e-10;
info=DSDPVMatShiftDiagonal(X,xshift); DSDPCHKBLOCKERR(blockj,info);
info=DSDPVMatScaleDiagonal(X,1.0+xscale); DSDPCHKBLOCKERR(blockj,info);
DSDPLogInfo(0,10,"XMat: shift diagonal: %4.2e, scale diagonal: %4.2e.\n",xshift,1+xscale);
}
DSDPFunctionReturn(0);
}
#undef __FUNCT__
#define __FUNCT__ "SDPConeComputeX"
int SDPConeComputeX(SDPCone sdpcone, int blockj, int n, double x[], int nn){
int info;
double mu=sdpcone->xmakermu;
double xnorm,xtrace,trxs;
char UPLQ;
DSDPVec DY=sdpcone->DYX,Y=sdpcone->YX,AX=sdpcone->Work;
DSDPVMat T;
DSDPFunctionBegin;
info=SDPConeCheckN(sdpcone,blockj,n);DSDPCHKBLOCKERR(blockj,info);
if (n<1){DSDPFunctionReturn(0);}
info=SDPConeGetStorageFormat(sdpcone,blockj,&UPLQ);DSDPCHKBLOCKERR(blockj,info);
info=DSDPMakeVMatWithArray(UPLQ,x,nn,n,&T);DSDPCHKBLOCKERR(blockj,info);
info=SDPConeComputeX3(sdpcone,blockj,mu,Y,DY,T);DSDPCHKBLOCKERR(blockj,info);
info=SDPConeComputeXDot(sdpcone,blockj,Y,T,AX,&xtrace,&xnorm,&trxs);DSDPCHKBLOCKERR(blockj,info);
info=DSDPVMatDestroy(&T);DSDPCHKBLOCKERR(blockj,info);
DSDPFunctionReturn(0);
}
#undef __FUNCT__
#define __FUNCT__ "SDPConeViewX"
int SDPConeViewX(SDPCone sdpcone, int blockj, int n, double x[], int nn){
int info;
char UPLQ;
DSDPVMat T;
DSDPFunctionBegin;
info=SDPConeCheckN(sdpcone,blockj,n);DSDPCHKBLOCKERR(blockj,info);
info=SDPConeGetStorageFormat(sdpcone,blockj,&UPLQ);DSDPCHKBLOCKERR(blockj,info);
info=DSDPMakeVMatWithArray(UPLQ,x,nn,n,&T);DSDPCHKBLOCKERR(blockj,info);
info=DSDPVMatView(T);DSDPCHKBLOCKERR(blockj,info);
info=DSDPVMatDestroy(&T);DSDPCHKBLOCKERR(blockj,info);
DSDPFunctionReturn(0);
}
#undef __FUNCT__
#define __FUNCT__ "SDPConeXVMultiply"
int SDPConeXVMultiply(SDPCone sdpcone, int blockj, double vin[], double vout[], int n){
int info;
SDPblk *blk=sdpcone->blk;
SDPConeVec V1,V2,V3,V4;
DSDPDualMat S,SS;
DSDPFunctionBegin;
info=SDPConeCheckN(sdpcone,blockj,n);DSDPCHKBLOCKERR(blockj,info);
if (sdpcone->blk[blockj].n>1){
S=blk[blockj].S; SS=blk[blockj].SS;
V2=blk[blockj].W; V3=blk[blockj].W2;
info=SDPConeVecCreateWArray(&V1,vin,n);
info=SDPConeVecCreateWArray(&V4,vout,n);
if (0){
info=DSDPDualMatCholeskySolveForward(S,V1,V3);DSDPCHKERR(info);
info=SDPConeVecScale(sqrt(sdpcone->xmakermu),V3);DSDPCHKERR(info);
info=DSDPDualMatCholeskySolveBackward(S,V3,V2);DSDPCHKERR(info);
info=DSDPDualMatCholeskyBackwardMultiply(SS,V2,V1);DSDPCHKERR(info);
}
info=DSDPDualMatCholeskyForwardMultiply(SS,V1,V2);DSDPCHKERR(info);
info=DSDPDualMatCholeskySolveForward(S,V2,V3);DSDPCHKERR(info);
info=SDPConeVecScale(sqrt(sdpcone->xmakermu),V3);DSDPCHKERR(info);
info=DSDPDualMatCholeskySolveBackward(S,V3,V4);DSDPCHKERR(info);
}
DSDPFunctionReturn(0);
}
#undef __FUNCT__
#define __FUNCT__ "SDPConeAddXVAV"
int SDPConeAddXVAV(SDPCone sdpcone, int blockj, double vin[], int n, double sum[], int mm){
int info;
SDPblk *blk=sdpcone->blk;
SDPConeVec V2;
DSDPVec W,Wout;
DSDPFunctionBegin;
info=SDPConeCheckN(sdpcone,blockj,n);DSDPCHKBLOCKERR(blockj,info);
W=sdpcone->Work;
info=DSDPVecSet(1.0,sdpcone->Work);DSDPCHKBLOCKERR(blockj,info);
if (sdpcone->blk[blockj].n>1){
info=SDPConeVecCreateWArray(&V2,vin,n);DSDPCHKERR(info);
info=DSDPVecCreateWArray(&Wout,sum,mm);DSDPCHKERR(info);
info=DSDPBlockvAv(&blk[blockj].ADATA,1.0,sdpcone->Work,V2,Wout);DSDPCHKBLOCKERR(blockj,info);
}
DSDPFunctionReturn(0);
}
#undef __FUNCT__
#define __FUNCT__ "SDPConeComputeXV"
int SDPConeComputeXV(SDPCone sdpcone, int blockj, int *derror){
int info; double rr;
DSDPVec Y,DY,W;
SDPblk *blk=sdpcone->blk;
DSDPTruth psdefinite1=DSDP_FALSE,psdefinite2=DSDP_FALSE;
DSDPDualMat S,SS;
DSDPVMat T;
DSDPFunctionBegin;
*derror=0;
info=SDPConeCheckJ(sdpcone,blockj);DSDPCHKBLOCKERR(blockj,info);
if (sdpcone->blk[blockj].n>1){
Y=sdpcone->YX; DY=sdpcone->DYX; W=sdpcone->Work;
T=blk[blockj].T; S=blk[blockj].S; SS=blk[blockj].SS;
info=DSDPVecWAXPY(W,-1.0,DY,Y);DSDPCHKBLOCKERR(blockj,info);
while (psdefinite1==DSDP_FALSE){
info=DSDPVecGetR(W,&rr);
info=DSDPVecSetR(W,10*rr-1e-12);
info=SDPConeComputeSS(sdpcone,blockj,W,T);DSDPCHKBLOCKERR(blockj,info);
info=DSDPDualMatSetArray(SS,T); DSDPCHKBLOCKERR(blockj,info);
info=DSDPDualMatCholeskyFactor(SS,&psdefinite1); DSDPCHKBLOCKERR(blockj,info);
}
while (psdefinite2==DSDP_FALSE){
info=SDPConeComputeSS(sdpcone,blockj,Y,T);DSDPCHKBLOCKERR(blockj,info);
info=DSDPDualMatSetArray(S,T); DSDPCHKBLOCKERR(blockj,info);
info=DSDPDualMatCholeskyFactor(S,&psdefinite2); DSDPCHKBLOCKERR(blockj,info);
if (psdefinite2==DSDP_FALSE){
info=DSDPVecGetR(Y,&rr);
info=DSDPVecSetR(Y,10*rr-1e-15);
}
}
if (psdefinite1==DSDP_FALSE || psdefinite2==DSDP_FALSE) *derror=1;
}
DSDPFunctionReturn(0);
}
SDPConeCheckJ
int SDPConeCheckJ(SDPCone sdpcone, int blockj)
Check validity of parameter.
Definition dsdpadddata.c:31
SDPConeCheckN
int SDPConeCheckN(SDPCone sdpcone, int blockj, int n)
Check validity of parameter.
Definition dsdpadddata.c:48
SDPConeCheckM
int SDPConeCheckM(SDPCone sdpcone, int m)
Check validity of parameter.
Definition dsdpadddata.c:68
DSDPTruth
DSDPTruth
Boolean variables.
Definition dsdpbasictypes.h:19
DSDP_FALSE
@ DSDP_FALSE
Definition dsdpbasictypes.h:19
DSDP_TRUE
@ DSDP_TRUE
Definition dsdpbasictypes.h:19
DSDPBlockASum
int DSDPBlockASum(DSDPBlockData *ADATA, double aa, DSDPVec Yk, DSDPVMat XX)
Sum the data matrices.
Definition dsdpblock.c:20
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
DSDPDualMatCholeskyFactor
int DSDPDualMatCholeskyFactor(DSDPDualMat S, DSDPTruth *psdefinite)
Factor the matrix.
Definition dsdpdualmat.c:320
DSDPDualMatCholeskyBackwardMultiply
int DSDPDualMatCholeskyBackwardMultiply(DSDPDualMat S, SDPConeVec B, SDPConeVec X)
Multiply by triangular matrix.
Definition dsdpdualmat.c:373
DSDPDualMatCholeskyForwardMultiply
int DSDPDualMatCholeskyForwardMultiply(DSDPDualMat S, SDPConeVec B, SDPConeVec X)
Multiply by triangular matrix.
Definition dsdpdualmat.c:346
DSDPDualMatCholeskySolveBackward
int DSDPDualMatCholeskySolveBackward(DSDPDualMat S, SDPConeVec B, SDPConeVec X)
Backward triangular solve.
Definition dsdpdualmat.c:295
DSDPDualMatCholeskySolveForward
int DSDPDualMatCholeskySolveForward(DSDPDualMat S, SDPConeVec B, SDPConeVec X)
Forward triangular solve.
Definition dsdpdualmat.c:267
DSDPDualMatIsFull
int DSDPDualMatIsFull(DSDPDualMat S, DSDPTruth *full)
Factor the matrix.
Definition dsdpdualmat.c:397
DSDPDualMatSetArray
int DSDPDualMatSetArray(DSDPDualMat S, DSDPVMat T)
Print the matrix.
Definition dsdpdualmat.c:160
DSDPDualMatInvert
int DSDPDualMatInvert(DSDPDualMat S)
Invert the matrix.
Definition dsdpdualmat.c:186
dsdpsdp.h
Internal SDPCone data structures and routines.
DSDPMakeVMatWithArray
int DSDPMakeVMatWithArray(char, double[], int, int, DSDPVMat *)
Allocate V matrix using the given array.
Definition sdpsss.c:381
SDPConeComputeXX
int SDPConeComputeXX(SDPCone, int, DSDPVec, double, DSDPDualMat, DSDPVMat)
Compute X.
Definition sdpcompute.c:235
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
DSDPVMatNormF2
int DSDPVMatNormF2(DSDPVMat X, double *normf2)
Compute square of Frobenius norm of matrix.
Definition dsdpxmat.c:186
DSDPVMatShiftDiagonal
int DSDPVMatShiftDiagonal(DSDPVMat X, double dadd)
Add something to diagonal elements.
Definition dsdpxmat.c:166
DSDPVMatZeroEntries
int DSDPVMatZeroEntries(DSDPVMat X)
Zero matrix.
Definition dsdpxmat.c:125
DSDPVMatDestroy
int DSDPVMatDestroy(DSDPVMat *X)
Deallocate matrix.
Definition dsdpxmat.c:86
DSDPVMatScaleDiagonal
int DSDPVMatScaleDiagonal(DSDPVMat X, double dscale)
Scaling diagonal is useful for inner products and norms.
Definition dsdpxmat.c:147
DSDPVMatView
int DSDPVMatView(DSDPVMat X)
Print matrix.
Definition dsdpxmat.c:107
SDPConeComputeX3
int SDPConeComputeX3(SDPCone sdpcone, int blockj, double mu, DSDPVec Y, DSDPVec DY, DSDPVMat X)
Compute the matrix X with the given information.
Definition sdpcone.c:140
SDPConeComputeSS
int SDPConeComputeSS(SDPCone sdpcone, int blockj, DSDPVec Y, DSDPVMat SS)
Sum the data matrices.
Definition sdpcone.c:18
SDPConeComputeXDot
int SDPConeComputeXDot(SDPCone sdpcone, int blockj, DSDPVec Y, DSDPVMat X, DSDPVec AX, double *xtrace, double *xnorm, double *tracexs)
Compute inner product of X with the Data, S, and norm of X.
Definition sdpcone.c:111
SDPConeVecScale
int SDPConeVecScale(double alpha, SDPConeVec VV)
Compute the Euclidean norm.
Definition sdpconevec.c:161
DSDPDualMat_C
Represents an S matrix for one block in the semidefinite cone.
Definition dsdpdualmat.h:18
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