00001 #include "dsdpsys.h"
00002 #include "dsdpdsmat_impl.h"
00003
00008 typedef struct {
00009 int n;
00010 double *an;
00011 int *col;
00012 int *nnz;
00013 } spdsmat;
00014
00015 static int SpSymMatSetURValuesP(void*DS, double v[], int nn, int n){
00016 spdsmat*ds=(spdsmat*)DS;
00017 int i,j,k1,k2,*nnz=ds->nnz,*col=ds->col;
00018 double *an=ds->an;
00019 for (i=0;i<n;i++,nnz++){
00020 k1=*nnz; k2=*(nnz+1);
00021 for (j=k1;j<k2;j++,an++,col++){
00022 if ((*col)==i){ *an = v[*col]/2;}
00023 else { *an = v[*col]; }
00024 }
00025 v+=i+1;
00026 }
00027 return 0;
00028 }
00029
00030 static int SpSymMatSetURValuesU(void*DS, double v[], int nn, int n){
00031 spdsmat*ds=(spdsmat*)DS;
00032 int i,j,k1,k2,*nnz=ds->nnz,*col=ds->col;
00033 double *an=ds->an;
00034 for (i=0;i<n;i++,nnz++){
00035 k1=*nnz; k2=*(nnz+1);
00036 for (j=k1;j<k2;j++,an++,col++){
00037 if ((*col)==i){ *an = v[*col]/2;}
00038 else { *an = v[*col]; }
00039 }
00040 v+=n;
00041 }
00042 return 0;
00043 }
00044
00045 static int SpSymMatView(void *DS){
00046 spdsmat*ds=(spdsmat*)DS;
00047 int i,j,k1,k2,n=ds->n,*nnz=ds->nnz,*col=ds->col;
00048 double *an=ds->an;
00049 for (i=0;i<n;i++){
00050 k1=nnz[i]; k2=nnz[i+1];
00051 printf("Row %d: ",i);
00052 for (j=k1;j<k2;j++){
00053 if (col[j]==i){ printf("%d: %4.4f",col[j],2*an[j]); }
00054 else { printf("%d: %4.4f",col[j],an[j]);}
00055 }
00056 printf("\n");
00057 }
00058 return 0;
00059 }
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071 static int SpSymMatDestroy(void *DS){
00072 spdsmat*ds=(spdsmat*)DS;
00073 int info;
00074 DSDPFREE(&ds->nnz,&info);if (info) return 1;
00075 DSDPFREE(&ds->col,&info);if (info) return 1;
00076 DSDPFREE(&ds->an,&info);if (info) return 1;
00077 DSDPFREE(&ds,&info);if (info) return 1;
00078 return 0;
00079 }
00080
00081 static int SpSymMatGetSize(void *DS, int*n){
00082 spdsmat*ds=(spdsmat*)DS;
00083 *n=ds->n;
00084 return 0;
00085 }
00086
00087 static int SpSymMatZero(void*DS){
00088 spdsmat*ds=(spdsmat*)DS;
00089 int nn=ds->nnz[ds->n];
00090 double *an=ds->an;
00091 memset((void*)an,0,nn*sizeof(double));
00092 return 0;
00093 }
00094
00095 static int SpSymMatMult(void*DS, double x[], double y[], int n){
00096 spdsmat*ds=(spdsmat*)DS;
00097 int i,j,k1,k2,*nnz=ds->nnz,*col=ds->col;
00098 double *an=ds->an;
00099 memset((void*)y,0,n*sizeof(double));
00100 for (i=0;i<n;i++,nnz++){
00101 k1=*nnz; k2=*(nnz+1);
00102 for (j=k1;j<k2;j++,col++,an++){
00103 y[*col] += x[i] * (*an);
00104 y[i] += x[*col] * (*an);
00105 }
00106 }
00107 return 0;
00108 }
00109
00110 static int SpSymMatVecVec(void*DS, double x[], int n, double *vAv){
00111 spdsmat*ds=(spdsmat*)DS;
00112 int i,j,k1,k2,*nnz=ds->nnz,*col=ds->col;
00113 double vv,*an=ds->an;
00114 *vAv=0;
00115 for (i=0;i<n;i++,nnz++){
00116 k1=*nnz; k2=*(nnz+1);
00117 vv=0;
00118 for (j=k1;j<k2;j++,col++,an++){
00119 vv+=x[*col]*(*an);
00120 }
00121 *vAv+=vv*x[i]*2;
00122 }
00123 return 0;
00124 }
00125
00126
00127
00128
00129
00130
00131
00132
00133
00134
00135
00136
00137
00138 static const char* dsmatname="SPARSE, SYMMETRIC MATRIX";
00139 static int DSDPDSSparseInitializeOpsP(struct DSDPDSMat_Ops* dsops){
00140 int info;
00141 if (!dsops) return 0;
00142 info=DSDPDSMatOpsInitialize(dsops); DSDPCHKERR(info);
00143 dsops->matseturmat=SpSymMatSetURValuesP;
00144 dsops->matview=SpSymMatView;
00145 dsops->matdestroy=SpSymMatDestroy;
00146 dsops->matgetsize=SpSymMatGetSize;
00147 dsops->matzeroentries=SpSymMatZero;
00148 dsops->matmult=SpSymMatMult;
00149 dsops->matvecvec=SpSymMatVecVec;
00150 dsops->id=6;
00151 dsops->matname=dsmatname;
00152 return 0;
00153 }
00154 static int DSDPDSSparseInitializeOpsU(struct DSDPDSMat_Ops* dsops){
00155 int info;
00156 if (!dsops) return 0;
00157 info=DSDPDSMatOpsInitialize(dsops); DSDPCHKERR(info);
00158 dsops->matseturmat=SpSymMatSetURValuesU;
00159 dsops->matview=SpSymMatView;
00160 dsops->matdestroy=SpSymMatDestroy;
00161 dsops->matgetsize=SpSymMatGetSize;
00162 dsops->matzeroentries=SpSymMatZero;
00163 dsops->matmult=SpSymMatMult;
00164 dsops->matvecvec=SpSymMatVecVec;
00165 dsops->id=6;
00166 dsops->matname=dsmatname;
00167 return 0;
00168 }
00169
00170 static struct DSDPDSMat_Ops tdsdsopsp;
00171 static struct DSDPDSMat_Ops tdsdsopsu;
00172 #undef __FUNCT__
00173 #define __FUNCT__ "DSDPCreateSparseDSMat"
00174 int DSDPSparseMatCreatePattern2P(int n, int rnnz[], int cols[], int tnnz,struct DSDPDSMat_Ops* *dsmatops, void**dsmat){
00175 int i,info;
00176 spdsmat*ds;
00177 DSDPFunctionBegin;
00178 DSDPCALLOC1(&ds,spdsmat,&info);DSDPCHKERR(info);
00179 DSDPCALLOC2(&ds->nnz,int,(n+1),&info);DSDPCHKERR(info);
00180 ds->nnz[0]=0;
00181 for (i=0;i<n;i++) ds->nnz[i+1]=ds->nnz[i]+rnnz[i];
00182 DSDPCALLOC2(&ds->col,int,tnnz,&info);DSDPCHKERR(info);
00183 DSDPCALLOC2(&ds->an,double,tnnz,&info);DSDPCHKERR(info);
00184 for (i=0;i<tnnz;i++) ds->col[i]=cols[i];
00185 info=DSDPDSSparseInitializeOpsP(&tdsdsopsp); DSDPCHKERR(info);
00186 *dsmatops=&tdsdsopsp;
00187 *dsmat=(void*)ds;
00188 DSDPFunctionReturn(0);
00189 }
00190
00191 #undef __FUNCT__
00192 #define __FUNCT__ "DSDPCreateSparseDSMatU"
00193 int DSDPSparseMatCreatePattern2U(int n, int rnnz[], int cols[], int tnnz,struct DSDPDSMat_Ops* *dsmatops, void**dsmat){
00194 int i,info;
00195 spdsmat*ds;
00196 DSDPFunctionBegin;
00197 DSDPCALLOC1(&ds,spdsmat,&info);DSDPCHKERR(info);
00198 DSDPCALLOC2(&ds->nnz,int,(n+1),&info);DSDPCHKERR(info);
00199 ds->nnz[0]=0;
00200 for (i=0;i<n;i++) ds->nnz[i+1]=ds->nnz[i]+rnnz[i];
00201 DSDPCALLOC2(&ds->col,int,tnnz,&info);DSDPCHKERR(info);
00202 DSDPCALLOC2(&ds->an,double,tnnz,&info);DSDPCHKERR(info);
00203 for (i=0;i<tnnz;i++) ds->col[i]=cols[i];
00204 info=DSDPDSSparseInitializeOpsU(&tdsdsopsu); DSDPCHKERR(info);
00205 *dsmatops=&tdsdsopsu;
00206 *dsmat=(void*)ds;
00207 DSDPFunctionReturn(0);
00208 }