spds.c

Go to the documentation of this file.
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 static int SpSymMatShiftDiagonal(void *DS, double dd){
00062   spdsmat*ds=(spdsmat*)DS;
00063   int i,n=ds->n,*nnz=ds->nnz;
00064   double *an=ds->an;
00065   for (i=0;i<n;i++){
00066     an[nnz[i+1]-1] += dd/2;
00067   }
00068   return 0;
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 static int SpSymMatAddRow(void *DS, int row, double dd, double v[], int n){
00127   spdsmat*ds=(spdsmat*)DS;
00128   int j,k1,k2,*nnz=ds->nnz,*col=ds->col;
00129   double *an=ds->an;
00130   k1=nnz[row]; k2=nnz[row+1];
00131   for (j=k1;j<k2;j++){
00132     if (row==col[j]){ an[j] += dd*v[col[j]]/2; } 
00133     else { an[j] += dd*v[col[j]]; }
00134   }
00135   return 0;
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 }

Generated on Wed May 21 19:31:19 2008 for DSDP by  doxygen 1.5.5