DSDP
dsdpdualmat.c
Go to the documentation of this file.
1 #include "dsdpdualmat_impl.h"
2 #include "dsdpdualmat.h"
3 #include "dsdpsys.h"
4 
10 #define DSDPNoOperationError(a); { DSDPSETERR1(1,"Dual natrix type: %s, Operation not defined\n",(a).dsdpops->matname);}
11 #define DSDPChkDMatError(a,b); { if (b){ DSDPSETERR1(b,"Dual natrix type: %s,\n",(a).dsdpops->matname);} }
12 
13 static int sdpdualsolve=0,sdpdualinvert=0;
14 
15 #undef __FUNCT__
16 #define __FUNCT__ "DSDPDualMatEventZero"
17 int DSDPDualMatEventZero(void){
18  DSDPFunctionBegin;
19  sdpdualinvert=0;sdpdualsolve=0;
20  DSDPFunctionReturn(0);
21 }
22 
23 #undef __FUNCT__
24 #define __FUNCT__ "DSDPDualMatEventInitialize"
25 int DSDPDualMatEventInitialize(void){
26  DSDPFunctionBegin;
27  if (sdpdualsolve==0){DSDPEventLogRegister("SDP SSolve",&sdpdualsolve);}
28  DSDPFunctionReturn(0);
29 }
30 
31 #undef __FUNCT__
32 #define __FUNCT__ "DSDPDualMatGetType"
33 int DSDPDualMatGetType(DSDPDualMat S, int *id){
34  DSDPFunctionBegin;
35  *id=S.dsdpops->id;
36  DSDPFunctionReturn(0);
37 }
38 
39 #undef __FUNCT__
40 #define __FUNCT__ "DSDPDualMatSetData"
49 int DSDPDualMatSetData(DSDPDualMat *S, struct DSDPDualMat_Ops* ops, void*data){
50  int info;
51  DSDPFunctionBegin;
52  (*S).dsdpops=ops;
53  (*S).matdata=data;
54  info=DSDPDualMatTest(*S);DSDPCHKERR(info);
55  DSDPFunctionReturn(0);
56 }
57 
58 #undef __FUNCT__
59 #define __FUNCT__ "DSDPDualMatDestroy"
66  int info;
67  DSDPFunctionBegin;
68  if ( S && (*S).dsdpops && (*S).dsdpops->matdestroy){
69  info=((*S).dsdpops->matdestroy)((*S).matdata); DSDPChkDMatError(*S,info);
70  } else {
71  /*
72  DSDPNoOperationError(*S);
73  */
74  }
75  info=DSDPDualMatSetData(S,0,0); DSDPCHKERR(info);
76  DSDPFunctionReturn(0);
77 }
78 
79 #undef __FUNCT__
80 #define __FUNCT__ "DSDPDualMatGetSize"
88  int info;
89  DSDPFunctionBegin;
90  if (S.dsdpops->matgetsize){
91  info=(S.dsdpops->matgetsize)(S.matdata,n); DSDPChkDMatError(S,info);
92  } else {
93  DSDPNoOperationError(S);
94  }
95  DSDPFunctionReturn(0);
96 }
97 
98 #undef __FUNCT__
99 #define __FUNCT__ "DSDPDualMatGetArray"
100 int DSDPDualMatGetArray(DSDPDualMat S, double **v, int *n){
101  int info;
102  DSDPFunctionBegin;
103  if (S.dsdpops->matgetarray){
104  info=(S.dsdpops->matgetarray)(S.matdata,v,n); DSDPChkDMatError(S,info);
105  } else {
106  *v=0;
107  *n=0;
108  }
109  DSDPFunctionReturn(0);
110 }
111 
112 
113 #undef __FUNCT__
114 #define __FUNCT__ "DSDPDualMatLogDeterminant"
123  int info;
124  DSDPFunctionBegin;
125  if (S.dsdpops->matlogdet){
126  info=(S.dsdpops->matlogdet)(S.matdata,logdet); DSDPChkDMatError(S,info);
127  } else {
128  DSDPNoOperationError(S);
129  }
130  DSDPFunctionReturn(0);
131 }
132 
133 #undef __FUNCT__
134 #define __FUNCT__ "DSDPDualMatView"
141  int info;
142  DSDPFunctionBegin;
143  if (S.dsdpops->matview){
144  info=(S.dsdpops->matview)(S.matdata); DSDPChkDMatError(S,info);
145  } else {
146  DSDPNoOperationError(S);
147  }
148  DSDPFunctionReturn(0);
149 }
150 
151 
152 #undef __FUNCT__
153 #define __FUNCT__ "DSDPDualMatSetArray"
161  double *ss;
162  int info,n,nn;
163  DSDPFunctionBegin;
164  if (S.dsdpops->matseturmat){
165  info=DSDPVMatGetSize(T,&n); DSDPCHKERR(info);
166  info=DSDPVMatGetArray(T,&ss,&nn); DSDPCHKERR(info);
167  info=(S.dsdpops->matseturmat)(S.matdata,ss,nn,n); DSDPChkDMatError(S,info);
168  info=DSDPVMatRestoreArray(T,&ss,&nn); DSDPCHKERR(info);
169  } else {
170  DSDPNoOperationError(S);
171  }
172  DSDPFunctionReturn(0);
173 }
174 
175 
176 #undef __FUNCT__
177 #define __FUNCT__ "DSDPDualMatInvert"
187  int info;
188  DSDPFunctionBegin;
189  /* DSDPEventLogBegin(sdpdualinvert); */
190  if (S.dsdpops->matinvert){
191  info=(S.dsdpops->matinvert)(S.matdata); DSDPChkDMatError(S,info);
192  } else {
193  DSDPNoOperationError(S);
194  }
195  /* DSDPEventLogEnd(sdpdualinvert); */
196  DSDPFunctionReturn(0);
197 }
198 
199 #undef __FUNCT__
200 #define __FUNCT__ "DSDPDualMatInverseAdd"
210  int info,n,nn;
211  double *ss;
212  DSDPFunctionBegin;
213  if (S.dsdpops->matinverseadd){
214  info=DSDPVMatGetSize(T,&n); DSDPCHKERR(info);
215  info=DSDPVMatGetArray(T,&ss,&nn); DSDPCHKERR(info);
216  info=(S.dsdpops->matinverseadd)(S.matdata,alpha,ss,nn,n); DSDPChkDMatError(S,info);
217  info=DSDPVMatRestoreArray(T,&ss,&nn); DSDPCHKERR(info);
218  } else {
219  DSDPNoOperationError(S);
220  }
221  DSDPFunctionReturn(0);
222 }
223 
224 
225 #undef __FUNCT__
226 #define __FUNCT__ "DSDPDualMatInverseMultiply"
237  int info,n;
238  double *bb,*xx;
239  DSDPFunctionBegin;
240  DSDPEventLogBegin(sdpdualsolve);
241  if (S.dsdpops->matinversemultiply){
242  info=SDPConeVecGetSize(X,&n); DSDPCHKERR(info);
243  info=SDPConeVecGetArray(B,&bb); DSDPCHKERR(info);
244  info=SDPConeVecGetArray(X,&xx); DSDPCHKERR(info);
245  info=(S.dsdpops->matinversemultiply)(S.matdata,IS.indx+1,IS.indx[0],bb,xx,n); DSDPChkDMatError(S,info);
246  info=SDPConeVecRestoreArray(X,&xx); DSDPCHKERR(info);
247  info=SDPConeVecRestoreArray(B,&bb); DSDPCHKERR(info);
248  } else {
249  DSDPNoOperationError(S);
250  }
251  DSDPEventLogEnd(sdpdualsolve);
252  DSDPFunctionReturn(0);
253 }
254 
255 
256 #undef __FUNCT__
257 #define __FUNCT__ "DSDPDualMatCholeskySolveForward"
268  int info,n;
269  double *bb,*xx;
270  DSDPFunctionBegin;
271  if (S.dsdpops->matsolveforward){
272  info=SDPConeVecGetSize(X,&n); DSDPCHKERR(info);
273  info=SDPConeVecGetArray(B,&bb); DSDPCHKERR(info);
274  info=SDPConeVecGetArray(X,&xx); DSDPCHKERR(info);
275  info=(S.dsdpops->matsolveforward)(S.matdata,bb,xx,n); DSDPChkDMatError(S,info);
276  info=SDPConeVecRestoreArray(X,&xx); DSDPCHKERR(info);
277  info=SDPConeVecRestoreArray(B,&bb); DSDPCHKERR(info);
278  } else {
279  DSDPNoOperationError(S);
280  }
281  DSDPFunctionReturn(0);
282 }
283 
284 #undef __FUNCT__
285 #define __FUNCT__ "DSDPDualMatDualMatCholeskySolveBackward"
296  int info,n;
297  double *bb,*xx;
298  DSDPFunctionBegin;
299  if (S.dsdpops->matsolvebackward){
300  info=SDPConeVecGetSize(X,&n); DSDPCHKERR(info);
301  info=SDPConeVecGetArray(B,&bb); DSDPCHKERR(info);
302  info=SDPConeVecGetArray(X,&xx); DSDPCHKERR(info);
303  info=(S.dsdpops->matsolvebackward)(S.matdata,bb,xx,n); DSDPChkDMatError(S,info);
304  info=SDPConeVecRestoreArray(X,&xx); DSDPCHKERR(info);
305  info=SDPConeVecRestoreArray(B,&bb); DSDPCHKERR(info);
306  } else {
307  DSDPNoOperationError(S);
308  }
309  DSDPFunctionReturn(0);
310 }
311 
312 #undef __FUNCT__
313 #define __FUNCT__ "DSDPDualMatCholeskyFactor"
321  int info;
322  int flag;
323  DSDPFunctionBegin;
324  if (S.dsdpops->matcholesky){
325  info=(S.dsdpops->matcholesky)(S.matdata,&flag); DSDPChkDMatError(S,info);
326  } else {
327  DSDPNoOperationError(S);
328  }
329  if (flag) *psdefinite=DSDP_FALSE;
330  else *psdefinite=DSDP_TRUE;
331  DSDPFunctionReturn(0);
332 }
333 
334 
335 #undef __FUNCT__
336 #define __FUNCT__ "DSDPDualMatCholeskyForwardMultiply"
347  int info,n;
348  double *bb,*xx;
349  DSDPFunctionBegin;
350  if (S.dsdpops->matforwardmultiply){
351  info=SDPConeVecGetSize(B,&n); DSDPCHKERR(info);
352  info=SDPConeVecGetArray(B,&bb); DSDPCHKERR(info);
353  info=SDPConeVecGetArray(X,&xx); DSDPCHKERR(info);
354  info=(S.dsdpops->matforwardmultiply)(S.matdata,bb,xx,n); DSDPChkDMatError(S,info);
355  info=SDPConeVecRestoreArray(X,&xx); DSDPCHKERR(info);
356  info=SDPConeVecRestoreArray(B,&bb); DSDPCHKERR(info);
357  } else {
358  DSDPNoOperationError(S);
359  }
360  DSDPFunctionReturn(0);
361 }
362 #undef __FUNCT__
363 #define __FUNCT__ "DSDPDualMatCholeskyBackwardMultiply"
374  int info,n;
375  double *bb,*xx;
376  DSDPFunctionBegin;
377  if (S.dsdpops->matbackwardmultiply){
378  info=SDPConeVecGetSize(B,&n); DSDPCHKERR(info);
379  info=SDPConeVecGetArray(B,&bb); DSDPCHKERR(info);
380  info=SDPConeVecGetArray(X,&xx); DSDPCHKERR(info);
381  info=(S.dsdpops->matbackwardmultiply)(S.matdata,bb,xx,n); DSDPChkDMatError(S,info);
382  info=SDPConeVecRestoreArray(X,&xx); DSDPCHKERR(info);
383  info=SDPConeVecRestoreArray(B,&bb); DSDPCHKERR(info);
384  } else {
385  DSDPNoOperationError(S);
386  }
387  DSDPFunctionReturn(0);
388 }
389 #undef __FUNCT__
390 #define __FUNCT__ "DSDPDualMatIsFull"
398  int info,flag=0;
399  DSDPFunctionBegin;
400  *full=DSDP_FALSE;
401  if (S.dsdpops->matfull){
402  info=(S.dsdpops->matfull)(S.matdata,&flag); DSDPChkDMatError(S,info);
403  } else {
404  DSDPNoOperationError(S);
405  }
406  if (flag) *full=DSDP_TRUE;
407  DSDPFunctionReturn(0);
408 }
409 
410 #undef __FUNCT__
411 #define __FUNCT__ "DSDPDataMatCheck"
412 int DSDPDualMatCheck(DSDPDualMat SS, SDPConeVec W1, SDPConeVec W2, DSDPIndex IS, DSDPVMat XX){
413  DSDPFunctionBegin;
414  DSDPFunctionReturn(0);
415 }
416 
417 static const char* dualmatname="NOT SET YET";
424  if (sops==NULL) return 0;
425  sops->matseturmat=0;
426  sops->matgetarray=0;
427  sops->matcholesky=0;
428  sops->matsolveforward=0;
429  sops->matsolvebackward=0;
430  sops->matinvert=0;
431  sops->matinverseadd=0;
432  sops->matinversemultiply=0;
433  sops->matforwardmultiply=0;
434  sops->matbackwardmultiply=0;
435  sops->matfull=0;
436  sops->matdestroy=0;
437  sops->matgetsize=0;
438  sops->matview=0;
439  sops->matlogdet=0;
440  sops->matname=dualmatname;
441  return 0;
442  }
443 
444 
445 static struct DSDPDualMat_Ops dsdpdualmatopsdefault;
446 
447 #undef __FUNCT__
448 #define __FUNCT__ "DSDPDualMatTest"
449 int DSDPDualMatTest(DSDPDualMat S){
450  int info;
451  DSDPFunctionBegin;
452  if (S.dsdpops==0 || S.dsdpops==&dsdpdualmatopsdefault){
453  } else if (S.dsdpops->mattest){
454  info=(S.dsdpops->mattest)(S.matdata); DSDPChkDMatError(S,info);
455  } else {
456  /*
457  DSDPNoOperationError(S);
458  */
459  }
460  DSDPFunctionReturn(0);
461 }
462 
463 
464 #undef __FUNCT__
465 #define __FUNCT__ "DSDPDualMatInitialize"
472  int info;
473  DSDPFunctionBegin;
474  info=DSDPDualMatOpsInitialize(&dsdpdualmatopsdefault);DSDPCHKERR(info);
475  info=DSDPDualMatSetData(S,&dsdpdualmatopsdefault,0); DSDPCHKERR(info);
476  DSDPFunctionReturn(0);
477 }
478 
DSDPTruth
Boolean variables.
@ DSDP_FALSE
@ DSDP_TRUE
int DSDPDualMatInitialize(DSDPDualMat *S)
Set pointers to null.
Definition: dsdpdualmat.c:471
int DSDPDualMatInverseAdd(DSDPDualMat S, double alpha, DSDPVMat T)
Add a multiple of the inverse to T.
Definition: dsdpdualmat.c:209
int DSDPDualMatView(DSDPDualMat S)
Print the matrix.
Definition: dsdpdualmat.c:140
int DSDPDualMatCholeskyFactor(DSDPDualMat S, DSDPTruth *psdefinite)
Factor the matrix.
Definition: dsdpdualmat.c:320
int DSDPDualMatCholeskyBackwardMultiply(DSDPDualMat S, SDPConeVec B, SDPConeVec X)
Multiply by triangular matrix.
Definition: dsdpdualmat.c:373
int DSDPDualMatCholeskyForwardMultiply(DSDPDualMat S, SDPConeVec B, SDPConeVec X)
Multiply by triangular matrix.
Definition: dsdpdualmat.c:346
int DSDPDualMatLogDeterminant(DSDPDualMat S, double *logdet)
Free the matrix structure.
Definition: dsdpdualmat.c:122
int DSDPDualMatOpsInitialize(struct DSDPDualMat_Ops *sops)
Set pointers to null.
Definition: dsdpdualmat.c:423
int DSDPDualMatCholeskySolveBackward(DSDPDualMat S, SDPConeVec B, SDPConeVec X)
Backward triangular solve.
Definition: dsdpdualmat.c:295
int DSDPDualMatDestroy(DSDPDualMat *S)
Free the matrix structure.
Definition: dsdpdualmat.c:65
int DSDPDualMatCholeskySolveForward(DSDPDualMat S, SDPConeVec B, SDPConeVec X)
Forward triangular solve.
Definition: dsdpdualmat.c:267
int DSDPDualMatIsFull(DSDPDualMat S, DSDPTruth *full)
Factor the matrix.
Definition: dsdpdualmat.c:397
int DSDPDualMatGetSize(DSDPDualMat S, int *n)
Free the matrix structure.
Definition: dsdpdualmat.c:87
int DSDPDualMatInverseMultiply(DSDPDualMat S, DSDPIndex IS, SDPConeVec B, SDPConeVec X)
Multiply the inverse by a vector or solve the system of equations.
Definition: dsdpdualmat.c:236
int DSDPDualMatSetArray(DSDPDualMat S, DSDPVMat T)
Print the matrix.
Definition: dsdpdualmat.c:160
int DSDPDualMatInvert(DSDPDualMat S)
Invert the matrix.
Definition: dsdpdualmat.c:186
The interface between the SDPCone and the matrix S.
Structure of function pointers that each symmetric positive definite matrix type (dense,...
Error handling, printing, and profiling.
int DSDPVMatGetArray(DSDPVMat X, double **v, int *nn)
Get the array that stores the matrix.
Definition: dsdpxmat.c:211
int DSDPVMatGetSize(DSDPVMat X, int *n)
Get number of rows and columns.
Definition: dsdpxmat.c:65
int DSDPVMatRestoreArray(DSDPVMat X, double **v, int *nn)
Restore the array that stores the matrix.
Definition: dsdpxmat.c:233
Represents an S matrix for one block in the semidefinite cone.
Definition: dsdpdualmat.h:18
Table of function pointers that operate on the S matrix.
Dense symmetric matrix for one block in the semidefinite cone.
Definition: dsdpxmat.h:17
Vector whose length corresponds to dimension of a block in a cone.
Definition: sdpconevec.h:13