DSDP
drowcol.c
Go to the documentation of this file.
1
2typedef struct {
3 int rc;
4 const double *val;
5 int n;
6 double x,y;
7} rcmat;
8
9#include "dsdpdatamat_impl.h"
10#include "dsdpsys.h"
17static int RCMatDestroy(void*);
18static int RCMatView(void*);
19static int RCMatVecVec(void*, double[], int, double *);
20static int RCMatDot(void*, double[], int, int, double *);
21static int RCMatGetRank(void*, int*, int);
22static int RCMatFactor(void*);
23static int RCMatGetEig(void*, int, double*, double[], int, int[], int*);
24static int RCMatAddRowMultiple(void*, int, double, double[], int);
25static int RCMatAddMultiple(void*, double, double[], int, int);
26static int RCMatGetRowNnz(void*, int, int[], int*, int);
27
28static int RCMatCreate(int n,int rowcol, const double v[], rcmat**M){
29 rcmat *M16;
30 M16=(rcmat*) malloc(1*sizeof(rcmat));
31 M16->val=v;
32 M16->rc=rowcol;
33 M16->n=n;
34 *M=M16;
35 return 0;
36}
37
38static int RCMatDestroy(void* AA){
39 free(AA);
40 return 0;
41}
42static int RCMatVecVec(void* A, double x[], int n, double *v){
43 rcmat*RC=(rcmat*)A;
44 int i;
45 double vv=0;
46 const double *val=RC->val;
47 for (i=0;i<n;i++){ vv+=val[i]*x[i];}
48 *v=2*vv*x[RC->rc];
49 return 0;
50}
51static int RCMatDot(void* A, double x[], int nn, int n1, double *v){
52 rcmat*RC=(rcmat*)A;
53 int i,k=0,rc=RC->rc,n=RC->n;
54 double vv=0;
55 const double *val=RC->val;
56 k=rc*(rc+1)/2;
57 for (i=0;i<=rc;i++,k++){ vv+=v[k]*val[i]; }
58 for (i=rc+1;i<n;i++,k+=i){ vv+=v[k+rc]*val[i];}
59 *v=2*vv;
60 return 0;
61}
62static int RCMatView(void* A){
63 rcmat*RC=(rcmat*)A;
64 int i;
65 printf("Row Col %d\n",RC->rc);
66 for (i=0;i<RC->n;i++){
67 printf("%4.4e ",RC->val[i]);
68 }
69 printf("\n");
70 return 0;
71}
72
73static int RCMatFactor(void* A){
74 /* Must compute x,y such that eigs are orthogonal */
75 rcmat*RC=(rcmat*)A;
76 int i,rc=RC->rc;
77 double vnorm2=0;
78 const double *val=RC->val;
79 for (i=0;i<RC->n;i++){ vnorm2+=val[i]*val[i];}
80 vnorm2=sqrt(vnorm2);
81 if (val[rc]>=0){
82 RC->y=-sqrt(vnorm2);
83 RC->x=sqrt(2*val[rc]+RC->y*RC->y);
84 } else {
85 RC->x=sqrt(vnorm2);
86 RC->y=-sqrt(-2*val[rc]+RC->x*RC->x);
87 }
88 return 0;
89}
90static int RCMatGetRank(void *A, int*rank, int n){
91 *rank=2;
92 return 0;
93}
94
95static int RCMatGetEig(void*A, int neig, double *eig, double v[], int n,int spind[], int *nind){
96 rcmat*RC=(rcmat*)A;
97 int i,rc=RC->rc;
98 double x=RC->x,y=RC->y,xmy=x-y;
99 const double *val=RC->val;
100 if (neig==0){
101 for (i=0;i<n;i++){spind[i]=i;}
102 for (i=0;i<n;i++){v[i]=val[i]/xmy;}
103 v[rc]=RC->x;
104 *nind=n;
105 *eig=1.0;
106 } else if (neig==1){
107 for (i=0;i<n;i++){spind[i]=i;}
108 for (i=0;i<n;i++){v[i]=val[i]/xmy;}
109 v[rc]=RC->y;
110 *nind=n;
111 *eig=-1.0;
112 } else { *eig=0;*nind=0;}
113 return 0;
114}
115
116static int RCMatGetRowNnz(void*A, int nrow, int nz[], int *nnzz, int n){
117 rcmat*RC=(rcmat*)A;
118 int i;
119 *nnzz=1;
120 if (nrow==RC->rc){for (i=0;i<n;i++){nz[i]++;}*nnzz=n;}
121 nz[nrow]++;
122 return 0;
123}
124
125static int RCMatAddRowMultiple(void*A, int nrow, double dd, double rrow[], int n){
126 rcmat*RC=(rcmat*)A;
127 int i;
128 if (nrow==RC->rc){
129 for (i=0;i<n;i++){rrow[i]+=dd*RC->val[i];}
130 }
131 rrow[nrow]+=dd*RC->val[nrow];
132 return 0;
133}
134static int RCMatAddMultiple(void*A, double dd, double vv[], int nn, int n1){
135 rcmat*RC=(rcmat*)A;
136 int i,rc=RC->rc,k=rc*(rc+1)/2,n=RC->n;
137 const double *val=RC->val;
138 for (i=0;i<=rc;i++,k++){
139 vv[k]+=dd*val[i];
140 }
141 for (i=rc+1;i<n;i++,k+=i){
142 vv[k+rc]+=dd*val[i];
143 }
144 return 0;
145}
146static int RCMatFNorm(void*A, int n, double *fnorm){
147 rcmat*RC=(rcmat*)A;
148 int i,rc=RC->rc;
149 double ff=0;
150 const double *val=RC->val;
151 for (i=0;i<n;i++){
152 ff+=val[i]*val[i];
153 }
154 ff*=2;
155 ff+=2*val[rc]*val[rc];
156 *fnorm=ff;
157 return 0;
158}
159static int RCMatCountNonzeros(void*A, int *nnz, int n){
160 *nnz=2*n-1;
161 return 0;
162}
163static struct DSDPDataMat_Ops rcmatops;
164static const char *datamatname="One Row and Column matrix";
165static int RCMatOperationsInitialize(struct DSDPDataMat_Ops* rcmatoperator){
166 int info;
167 if (rcmatoperator==NULL) return 0;
168 info=DSDPDataMatOpsInitialize(rcmatoperator); if (info){ return info;}
169 rcmatoperator->matfactor1=RCMatFactor;
170 rcmatoperator->matgetrank=RCMatGetRank;
171 rcmatoperator->matgeteig=RCMatGetEig;
172 rcmatoperator->matvecvec=RCMatVecVec;
173 rcmatoperator->matrownz=RCMatGetRowNnz;
174 rcmatoperator->matdot=RCMatDot;
175 rcmatoperator->matfnorm2=RCMatFNorm;
176 rcmatoperator->matnnz=RCMatCountNonzeros;
177 rcmatoperator->mataddrowmultiple=RCMatAddRowMultiple;
178 rcmatoperator->mataddallmultiple=RCMatAddMultiple;
179 rcmatoperator->matdestroy=RCMatDestroy;
180 rcmatoperator->matview=RCMatView;
181 rcmatoperator->matname=datamatname;
182 rcmatoperator->id=27;
183 return 0;
184}
185
186
187#undef __FUNCT__
188#define __FUNCT__ "DSDPGetRCMat"
189int DSDPGetRCMat(int n,const double val[],int rowcol, struct DSDPDataMat_Ops**sops, void**smat){
190 rcmat *AA;
191 int info;
192 DSDPFunctionBegin;
193 info=RCMatCreate(n,rowcol,val,&AA);DSDPCHKERR(info);
194 info=RCMatOperationsInitialize(&rcmatops); DSDPCHKERR(info);
195 if (sops){*sops=&rcmatops;}
196 if (smat){*smat=(void*)AA;}
197 DSDPFunctionReturn(0);
198}
199
int DSDPDataMatOpsInitialize(struct DSDPDataMat_Ops *dops)
Initialize the table of function pointers for SDP Data matrices.
Definition dsdpdatamat.c:47
Structure of function pointers that each SDP data matrix type (sparse, dense, constant,...
Error handling, printing, and profiling.
Table of function pointers that operate on the data matrix.