00001 #include "dsdpcone_impl.h"
00002 #include "dsdpsys.h"
00003 #include "dsdpbasictypes.h"
00010 struct RDCone{
00011 double primalr,dualr,x,logr;
00012 DSDPPenalty UsePenalty;
00013 DSDP dsdp;
00014 };
00015
00016 typedef struct RDCone RCone;
00017
00018 #undef __FUNCT__
00019 #define __FUNCT__ "DSDPRHessian"
00020 static int DSDPRHessian( void *dspcone, double mu, DSDPSchurMat M, DSDPVec vrhs1, DSDPVec vrhs2){
00021 RCone *K=(RCone*)dspcone;
00022 int info,m;
00023 double rr,sr,rssr;
00024 DSDPFunctionBegin;
00025 if (K->dualr ){
00026 info=DSDPVecGetSize(vrhs2,&m);DSDPCHKERR(info);
00027 info=DSDPSchurMatVariableCompute(M,m-1,&rr);DSDPCHKERR(info);
00028 if (rr){
00029 sr=-mu*rr/(K->dualr);
00030 rssr=mu*rr/(K->dualr*K->dualr);
00031 info=DSDPVecAddR(vrhs2,sr);DSDPCHKERR(info);
00032 info=DSDPSchurMatAddDiagonalElement(M,m-1,rssr);DSDPCHKERR(info);
00033 }
00034 }
00035 DSDPFunctionReturn(0);
00036 }
00037
00038
00039 #undef __FUNCT__
00040 #define __FUNCT__ "DSDPRHS"
00041 static int DSDPRHS( void *dspcone, double mu, DSDPVec vrow,DSDPVec vrhs1,DSDPVec vrhs2){
00042 RCone *K=(RCone*)dspcone;
00043 int info;
00044 double rr,sr;
00045 DSDPFunctionBegin;
00046 if (K->dualr ){
00047 sr=-mu/(K->dualr);
00048 info=DSDPVecGetR(vrow,&rr);DSDPCHKERR(info);
00049 info=DSDPVecAddR(vrhs2,rr*sr);DSDPCHKERR(info);
00050 }
00051 DSDPFunctionReturn(0);
00052 }
00053
00054
00055 #undef __FUNCT__
00056 #define __FUNCT__ "DSDPSetupRCone"
00057 static int DSDPSetupRCone(void* dspcone,DSDPVec y){
00058 DSDPFunctionBegin;
00059 DSDPFunctionReturn(0);
00060 }
00061
00062 #undef __FUNCT__
00063 #define __FUNCT__ "DSDPSetupRCone2"
00064 static int DSDPSetupRCone2(void* dspcone, DSDPVec y, DSDPSchurMat M){
00065 DSDPFunctionBegin;
00066 DSDPFunctionReturn(0);
00067 }
00068
00069
00070 #undef __FUNCT__
00071 #define __FUNCT__ "DSDPDestroyRCone"
00072 static int DSDPDestroyRCone(void* dspcone){
00073 int info;
00074 DSDPFunctionBegin;
00075 DSDPFREE(&dspcone,&info);DSDPCHKERR(info);
00076 DSDPFunctionReturn(0);
00077 }
00078
00079
00080 #undef __FUNCT__
00081 #define __FUNCT__ "DSDPRSize"
00082 static int DSDPRSize(void*dspcone,double *n){
00083 RCone *K=(RCone*)dspcone;
00084 DSDPFunctionBegin;
00085 if (K->dualr){*n=1;}
00086 else {*n=0;}
00087 DSDPFunctionReturn(0);
00088 }
00089
00090 #undef __FUNCT__
00091 #define __FUNCT__ "DSDPRSparsity"
00092 static int DSDPRSparsity(void*dspcone,int row, int *tnnz, int rnnz[], int m){
00093 DSDPFunctionBegin;
00094 *tnnz=0;
00095 DSDPFunctionReturn(0);
00096 }
00097
00098
00099 #undef __FUNCT__
00100 #define __FUNCT__ "DSDPComputeRS"
00101 static int DSDPComputeRS(void *dspcone, DSDPVec Y, DSDPDualFactorMatrix flag, DSDPTruth *ispsdefinite){
00102 RCone *K=(RCone*)dspcone;
00103 int info;
00104 double rbeta;
00105 DSDPFunctionBegin;
00106 info=DSDPVecGetR(Y,&rbeta); DSDPCHKERR(info);
00107 if (K->UsePenalty==DSDPAlways){
00108 if (rbeta<0){ *ispsdefinite=DSDP_TRUE; } else { *ispsdefinite=DSDP_FALSE;}
00109 } else {
00110 if (rbeta>0) rbeta=0;
00111 *ispsdefinite=DSDP_TRUE;
00112 }
00113 if (flag==DUAL_FACTOR){ K->dualr=rbeta; } else { K->primalr=rbeta;}
00114 DSDPFunctionReturn(0);
00115 }
00116 #undef __FUNCT__
00117 #define __FUNCT__ "DSDPInvertRS"
00118 static int DSDPInvertRS(void *dspcone){
00119 DSDPFunctionBegin;
00120 DSDPFunctionReturn(0);
00121 }
00122
00123
00124
00125 #undef __FUNCT__
00126 #define __FUNCT__ "DSDPComputeRStepLength"
00127 static int DSDPComputeRStepLength(void *dspcone, DSDPVec DY, DSDPDualFactorMatrix flag, double *maxsteplength){
00128 RCone *K=(RCone*)dspcone;
00129 double r,rbeta,msteplength=1.0e100,rt=1.0e30;
00130 int info;
00131 DSDPFunctionBegin;
00132
00133 info=DSDPVecGetR(DY,&rbeta); DSDPCHKERR(info);
00134 if (flag==DUAL_FACTOR){ r=K->dualr; } else { r=K->primalr;}
00135 if (r * rbeta<0) rt=-r/rbeta;
00136
00137 if (K->UsePenalty==DSDPAlways){msteplength=rt;}
00138 else if (flag==PRIMAL_FACTOR){ msteplength=rt;}
00139 else if (flag==DUAL_FACTOR){msteplength=rt/0.94;}
00140
00141 *maxsteplength=msteplength;
00142 DSDPFunctionReturn(0);
00143 }
00144
00145 #undef __FUNCT__
00146 #define __FUNCT__ "DSDPRX"
00147 static int DSDPRX( void *dspcone, double mu, DSDPVec y, DSDPVec dy, DSDPVec AX,double *tracexs){
00148 RCone *K=(RCone*)dspcone;
00149 int info;
00150 double rr,dr,trxs,r=K->dualr;
00151 DSDPFunctionBegin;
00152
00153 info=DSDPVecGetR(y,&rr); DSDPCHKERR(info);
00154 info=DSDPVecGetR(dy,&dr); DSDPCHKERR(info);
00155 if (K->dualr){
00156 r=-1.0/rr;
00157 K->x=mu*(r-r*dr*r);
00158 trxs=K->x/r;
00159 DSDPLogInfo(0,2,"RESIDUAL X (Minimum Penalty Parameter): %4.4e, Trace(XS): %4.4e\n",K->x,trxs);
00160
00161 } else {
00162 K->x=0.0;
00163 }
00164 DSDPFunctionReturn(0);
00165 }
00166
00167 #undef __FUNCT__
00168 #define __FUNCT__ "DSDPSetX"
00169 static int DSDPSetX( void *dspcone, double mu, DSDPVec y, DSDPVec dy){
00170 DSDPFunctionBegin;
00171
00172
00173
00174
00175
00176
00177
00178 DSDPFunctionReturn(0);
00179 }
00180
00181 #undef __FUNCT__
00182 #define __FUNCT__ "DSDPComputeRLog"
00183 static int DSDPComputeRLog(void *dspcone, double *logobj, double *logdet){
00184 RCone *K=(RCone*)dspcone;
00185 DSDPFunctionBegin;
00186 *logdet=K->logr;
00187 *logobj=0;
00188 if (K->dualr<0){
00189 *logdet=log(-K->dualr);
00190 K->logr=log(-K->dualr);
00191 }
00192 DSDPFunctionReturn(0);
00193 }
00194
00195 #undef __FUNCT__
00196 #define __FUNCT__ "DSDPRANorm2"
00197 static int DSDPRANorm2(void *dspcone,DSDPVec Anorm2){
00198 DSDPFunctionBegin;
00199 DSDPFunctionReturn(0);
00200 }
00201
00202 #undef __FUNCT__
00203 #define __FUNCT__ "DSDPRMultiplyAdd"
00204 static int DSDPRMultiplyAdd(void *dspcone,double mu,DSDPVec vrow,DSDPVec vin,DSDPVec vout){
00205 RCone *K=(RCone*)dspcone;
00206 int info;
00207 double v1,v2,rssr;
00208 DSDPFunctionBegin;
00209 if (K->dualr){
00210 info=DSDPVecGetR(vrow,&v1);DSDPCHKERR(info);
00211 info=DSDPVecGetR(vin,&v2);DSDPCHKERR(info);
00212 rssr=v1*v2*mu/(K->dualr*K->dualr);
00213 info=DSDPVecAddR(vout,rssr);DSDPCHKERR(info);
00214 }
00215 DSDPFunctionReturn(0);
00216 }
00217
00218 #undef __FUNCT__
00219 #define __FUNCT__ "DSDPRMonitor"
00220 static int DSDPRMonitor( void *dspcone, int tag){
00221 DSDPFunctionBegin;
00222 DSDPFunctionReturn(0);
00223 }
00224
00225 static struct DSDPCone_Ops kops;
00226 static const char* matname="R Cone";
00227
00228 #undef __FUNCT__
00229 #define __FUNCT__ "RConeOperationsInitialize"
00230 static int RConeOperationsInitialize(struct DSDPCone_Ops* coneops){
00231 int info;
00232 if (coneops==NULL) return 0;
00233 info=DSDPConeOpsInitialize(coneops); DSDPCHKERR(info);
00234 coneops->conehessian=DSDPRHessian;
00235 coneops->conesetup=DSDPSetupRCone;
00236 coneops->conesetup2=DSDPSetupRCone2;
00237 coneops->conedestroy=DSDPDestroyRCone;
00238 coneops->conecomputes=DSDPComputeRS;
00239 coneops->coneinverts=DSDPInvertRS;
00240 coneops->conesetxmaker=DSDPSetX;
00241 coneops->conecomputex=DSDPRX;
00242 coneops->conerhs=DSDPRHS;
00243 coneops->conemaxsteplength=DSDPComputeRStepLength;
00244 coneops->conelogpotential=DSDPComputeRLog;
00245 coneops->conesize=DSDPRSize;
00246 coneops->conesparsity=DSDPRSparsity;
00247 coneops->coneanorm2=DSDPRANorm2;
00248 coneops->conemonitor=DSDPRMonitor;
00249 coneops->conehmultiplyadd=DSDPRMultiplyAdd;
00250 coneops->id=19;
00251 coneops->name=matname;
00252 return 0;
00253 }
00254
00255
00256
00257
00258
00259
00260
00261
00262 #undef __FUNCT__
00263 #define __FUNCT__ "RConeSetType"
00264 int RConeSetType(RCone *rcone, DSDPPenalty UsePenalty){
00265 DSDPFunctionBegin;
00266 rcone->UsePenalty=UsePenalty;
00267 DSDPFunctionReturn(0);
00268 }
00269
00277 #undef __FUNCT__
00278 #define __FUNCT__ "DSDPAddRCone"
00279 int DSDPAddRCone(DSDP dsdp, RCone **rrcone){
00280 RCone *rcone;
00281 DSDPPenalty UsePenalty=DSDPInfeasible;
00282 int info;
00283 DSDPFunctionBegin;
00284 info=RConeOperationsInitialize(&kops); DSDPCHKERR(info);
00285 DSDPCALLOC1(&rcone,RCone,&info); DSDPCHKERR(info);
00286 info=RConeSetType(rcone,UsePenalty); DSDPCHKERR(info);
00287 rcone->dsdp=dsdp;
00288 rcone->logr=0;
00289 *rrcone=rcone;
00290 info=DSDPAddCone(dsdp,&kops,(void*)rcone); DSDPCHKERR(info);
00291 DSDPFunctionReturn(0);
00292 }
00293