DSDP
dsdpschurmat.c
Go to the documentation of this file.
1#include "dsdpschurmat_impl.h"
2#include "dsdpschurmat.h"
3#include "dsdpbasictypes.h"
4#include "dsdpsys.h"
5
11static int hfactorevent=0,hsolveevent=0;
12
13#define DSDPNoOperationError(a); { DSDPSETERR1(10,"Schur matrix type: %s, Operation not defined\n",(a).dsdpops->matname); }
14#define DSDPChkMatError(a,b); { if (b){ DSDPSETERR1(b,"Schur matrix type: %s,\n",(a).dsdpops->matname);} }
15
16static int DSDPApplySMW(DSDPSchurMat, DSDPVec, DSDPVec);
17static int DSDPSchurMatSolveM(DSDPSchurMat, DSDPVec, DSDPVec);
18
19#undef __FUNCT__
20#define __FUNCT__ "DSDPSchurMatSetData"
28int DSDPSchurMatSetData(DSDPSchurMat *M, struct DSDPSchurMat_Ops* ops, void*data){
29 DSDPFunctionBegin;
30 (*M).dsdpops=ops;
31 (*M).data=data;
32 DSDPFunctionReturn(0);
33}
34
35static const char* schurmatname="NOT NAMED YET";
36
37#undef __FUNCT__
38#define __FUNCT__ "DSDPSchurMatOpsInitialize"
44int DSDPSchurMatOpsInitialize(struct DSDPSchurMat_Ops* dops){
45 DSDPFunctionBegin;
46 if (dops==NULL) return 0;
47 dops->matzero=0;
48 dops->matrownonzeros=0;
49 dops->mataddrow=0;
50 dops->mataddelement=0;
51 dops->matadddiagonal=0;
52 dops->matshiftdiagonal=0;
53 dops->matassemble=0;
54 dops->matscaledmultiply=0;
55 dops->matmultr=0;
56 dops->matfactor=0;
57 dops->matsolve=0;
58 dops->pmatonprocessor=0;
59 dops->pmatwhichdiag=0;
60 dops->pmatdistributed=0;
61 dops->matdestroy=0;
62 dops->matview=0;
63 dops->matsetup=0;
64 dops->id=0;
65 dops->matname=schurmatname;
66 DSDPFunctionReturn(0);
67}
68
69static struct DSDPSchurMat_Ops dsdpmops;
70
71
72#undef __FUNCT__
73#define __FUNCT__ "DSDPSchurMatOpsInitialize"
80 int info;
81 DSDPFunctionBegin;
82 info=DSDPSchurMatOpsInitialize(&dsdpmops); DSDPCHKERR(info);
83 info=DSDPSchurMatSetData(M,&dsdpmops,0); DSDPCHKERR(info);
84 DSDPCALLOC1(&M->schur,DSDPSchurInfo,&info);DSDPCHKERR(info);
85 M->schur->m=0; M->schur->r=0; M->schur->dd=0;
86 info=DSDPInitializeFixedVariable(&M->schur->fv);DSDPCHKERR(info);
87 DSDPFunctionReturn(0);
88}
89
90#undef __FUNCT__
91#define __FUNCT__ "DSDPSchurMatZeroEntries"
98 int info;
99 DSDPFunctionBegin;
100 if (M.dsdpops->matzero){
101 info=(M.dsdpops->matzero)(M.data); DSDPChkMatError(M,info);
102 } else {
103 DSDPNoOperationError(M);
104 }
105 DSDPFunctionReturn(0);
106}
107
108
109#undef __FUNCT__
110#define __FUNCT__ "DSDPSchurMatShiftDiagonal"
121 int info;
122 DSDPFunctionBegin;
123 if (dd==0){DSDPFunctionReturn(0);}
124 M.schur->dd=dd;
125 if (M.dsdpops->matshiftdiagonal){
126 /* if(M.schur->r){info=DSDPVecAddR(M.schur->rhs3,dd);DSDPCHKERR(info);} */
127 info=(M.dsdpops->matshiftdiagonal)(M.data,dd); DSDPChkMatError(M,info);
128 DSDPLogInfo(0,2,"Add %4.4e to the Diagonal of Schur Matrix\n",dd);
129 } else {
130 DSDPNoOperationError(M);
131 }
132 DSDPFunctionReturn(0);
133}
134
135
136#undef __FUNCT__
137#define __FUNCT__ "DSDPSchurMatInParallel"
150 int info,flg;
151 DSDPFunctionBegin;
152 if (M.dsdpops->pmatdistributed){
153 info=(M.dsdpops->pmatdistributed)(M.data,&flg); DSDPChkMatError(M,info);
154 if (flg) *flag=DSDP_TRUE; else *flag=DSDP_FALSE;
155 } else {
156 *flag=DSDP_FALSE;
157 }
158 DSDPFunctionReturn(0);
159}
160
161
162
163#undef __FUNCT__
164#define __FUNCT__ "DSDPSchurMatAssemble"
175 int info;
176 DSDPFunctionBegin;
177 if (M.dsdpops->matassemble){
178 info=(M.dsdpops->matassemble)(M.data); DSDPChkMatError(M,info);
179 } else {
180 DSDPNoOperationError(M);
181 }
182 DSDPFunctionReturn(0);
183}
184
185#undef __FUNCT__
186#define __FUNCT__ "DSDPSchurMatFactor"
197 int info,flag=0;
198 DSDPVec rhs3=M.schur->rhs3,dy3=M.schur->dy3;
199 DSDPFunctionBegin;
200 *successful=DSDP_TRUE;
201 DSDPEventLogBegin(hfactorevent);
202 if (M.dsdpops->matfactor){
203 info=(M.dsdpops->matfactor)(M.data,&flag); DSDPChkMatError(M,info);
204 if (flag){
205 *successful=DSDP_FALSE;
206 DSDPLogInfo(0,2,"Indefinite Schur Matrix -- Bad Factorization\n");
207 }
208 } else {
209 DSDPNoOperationError(M);
210 }
211 DSDPEventLogEnd(hfactorevent);
212 if (M.schur->r){
213 info=DSDPSchurMatSolveM(M,rhs3,dy3);DSDPCHKERR(info);}
214 else {info=DSDPVecZero(dy3);DSDPCHKERR(info);}
215 DSDPFunctionReturn(0);
216}
217
218
219#undef __FUNCT__
220#define __FUNCT__ "DSDPSchurMatMultiply"
232 int info,n;
233 double *xx,*yy,r=M.schur->r;
234 double r1,r2,dd;
235 DSDPVec rhs3;
236 DSDPFunctionBegin;
237
238 if (M.dsdpops->matscaledmultiply){
239 info=DSDPVecGetSize(x,&n); DSDPCHKERR(info);
240 info=DSDPVecGetArray(x,&xx); DSDPCHKERR(info);
241 info=DSDPVecGetArray(y,&yy); DSDPCHKERR(info);
242 info=(M.dsdpops->matscaledmultiply)(M.data,xx+1,yy+1,n-2); DSDPChkMatError(M,info);
243 yy[0]=0;
244 yy[n-1]=0;
245 info=DSDPVecRestoreArray(y,&yy); DSDPCHKERR(info);
246 info=DSDPVecRestoreArray(x,&xx); DSDPCHKERR(info);
247 } else {
248 DSDPNoOperationError(M);
249 }
250 if (r){
251 rhs3=M.schur->rhs3;
252 info=DSDPVecGetR(rhs3,&r2);DSDPCHKERR(info);
253 info=DSDPVecGetR(x,&r1);DSDPCHKERR(info);
254 info=DSDPVecAXPY(r1,rhs3,y);DSDPCHKERR(info);
255 info=DSDPVecDot(rhs3,x,&dd);DSDPCHKERR(info);
256 info=DSDPVecAddR(y,dd-r1*r2);DSDPCHKERR(info);
257 }
258 DSDPFunctionReturn(0);
259}
260
261#undef __FUNCT__
262#define __FUNCT__ "DSDPSchurMatMultR"
263int DSDPSchurMatMultR(DSDPSchurMat M, DSDPVec x, DSDPVec y){
264 int info,n;
265 double *xx,*yy,r=M.schur->r;
266 double r1,r2,dd;
267 DSDPVec rhs3;
268 DSDPFunctionBegin;
269
270 if (M.dsdpops->matmultr){
271 info=DSDPVecGetSize(x,&n); DSDPCHKERR(info);
272 info=DSDPVecGetArray(x,&xx); DSDPCHKERR(info);
273 info=DSDPVecGetArray(y,&yy); DSDPCHKERR(info);
274 info=(M.dsdpops->matmultr)(M.data,xx+1,yy+1,n-2); DSDPChkMatError(M,info);
275 yy[0]=0;
276 yy[n-1]=0;
277 info=DSDPVecRestoreArray(y,&yy); DSDPCHKERR(info);
278 info=DSDPVecRestoreArray(x,&xx); DSDPCHKERR(info);
279 if (r){
280 rhs3=M.schur->rhs3;
281 info=DSDPVecGetR(rhs3,&r2);DSDPCHKERR(info);
282 info=DSDPVecGetR(x,&r1);DSDPCHKERR(info);
283 info=DSDPVecAXPY(r1,rhs3,y);DSDPCHKERR(info);
284 info=DSDPVecDot(rhs3,x,&dd);DSDPCHKERR(info);
285 info=DSDPVecAddR(y,dd-r1*r2);DSDPCHKERR(info);
286 }
287 } else {
288 info=DSDPVecZero(y);DSDPCHKERR(info);
289 /* DSDPNoOperationError(M); */
290 }
291 DSDPFunctionReturn(0);
292}
293
294#undef __FUNCT__
295#define __FUNCT__ "DSDPSchurMatReducePVec"
308 int info,n;
309 double *xx;
310 DSDPTruth flag;
311 DSDPFunctionBegin;
312
313 if (M.dsdpops->pmatreduction){
314 info=DSDPVecGetSize(x,&n); DSDPCHKERR(info);
315 info=DSDPVecGetArray(x,&xx); DSDPCHKERR(info);
316 info=(M.dsdpops->pmatreduction)(M.data,xx+1,n-2); DSDPChkMatError(M,info);
317 info=DSDPVecRestoreArray(x,&xx); DSDPCHKERR(info);
318 } else {
319 info=DSDPSchurMatInParallel(M,&flag);DSDPChkMatError(M,info);
320 if (flag==DSDP_TRUE){
321 DSDPNoOperationError(M);
322 }
323 }
324 info=DSDPZeroFixedVariables(M,x);DSDPCHKERR(info);
325 DSDPFunctionReturn(0);
326}
327
328
329
330#undef __FUNCT__
331#define __FUNCT__ "DSDPSchurMatSetR"
339 DSDPFunctionBegin;
340 M.schur->r=rr;
341 DSDPFunctionReturn(0);
342}
343
344#undef __FUNCT__
345#define __FUNCT__ "DSDPSchurMatSetup"
353 int info,m;
354 DSDPFunctionBegin;
355 info=DSDPVecDuplicate(Y,&M.schur->rhs3);
356 info=DSDPVecDuplicate(Y,&M.schur->dy3);
357 info=DSDPVecGetSize(Y,&m);DSDPCHKERR(info);
358 if (M.dsdpops->matsetup){
359 info=(M.dsdpops->matsetup)(M.data,m-2); DSDPChkMatError(M,info);
360 } else {
361 DSDPNoOperationError(M);
362 }
363 DSDPEventLogRegister("Factor Newton Eq.",&hfactorevent);
364 DSDPEventLogRegister("Solve Newton Eq.",&hsolveevent);
365 DSDPFunctionReturn(0);
366}
367
368
369#undef __FUNCT__
370#define __FUNCT__ "DSDPSchurMatView"
377 int info;
378 DSDPFunctionBegin;
379 if (M.dsdpops->matview){
380 info=(M.dsdpops->matview)(M.data); DSDPChkMatError(M,info);
381 } else {
382 DSDPNoOperationError(M);
383 }
384 info=DSDPVecView(M.schur->rhs3);DSDPCHKERR(info);
385 DSDPFunctionReturn(0);
386}
387
388
389#undef __FUNCT__
390#define __FUNCT__ "DSDPSchurMatRowScaling"
400 int info;
401 DSDPFunctionBegin;
402 info=DSDPSchurMatDiagonalScaling(M,D);DSDPCHKERR(info);
403 info=DSDPZeroFixedVariables(M,D);DSDPCHKERR(info);
404 DSDPFunctionReturn(0);
405}
406
407#undef __FUNCT__
408#define __FUNCT__ "DSDPSchurMatDestroy"
415 int info;
416 DSDPFunctionBegin;
417 if ((*M).dsdpops->matdestroy){
418 info=((*M).dsdpops->matdestroy)((*M).data); DSDPChkMatError(*M,info);
419 } else {
420 /*
421 DSDPNoOperationError(*M);
422 */
423 }
424 info=DSDPVecDestroy(&M->schur->rhs3);DSDPCHKERR(info);
425 info=DSDPVecDestroy(&M->schur->dy3);DSDPCHKERR(info);
426 /* info=DSDPSchurMatSetData(M,0,0); DSDPCHKERR(info); */
427 info=DSDPSchurMatOpsInitialize(&dsdpmops); DSDPCHKERR(info);
428 info=DSDPSchurMatSetData(M,&dsdpmops,0); DSDPCHKERR(info);
429 DSDPFREE(&M->schur,&info);DSDPCHKERR(info);
430 DSDPFunctionReturn(0);
431}
432
433#undef __FUNCT__
434#define __FUNCT__ "DSDPSchurMatSolveM"
435static int DSDPSchurMatSolveM(DSDPSchurMat M, DSDPVec b, DSDPVec x){
436 int info,n;
437 double *xx,*bb;
438 DSDPFunctionBegin;
439 info=DSDPEventLogBegin(hsolveevent);
440 if (M.dsdpops->matsolve){
441 info=DSDPVecGetArray(b,&bb); DSDPCHKERR(info);
442 info=DSDPVecGetSize(x,&n); DSDPCHKERR(info);
443 info=DSDPVecZero(x);DSDPCHKERR(info);
444 info=DSDPVecGetArray(x,&xx); DSDPCHKERR(info);
445 info=(M.dsdpops->matsolve)(M.data,bb+1,xx+1,n-2); DSDPChkMatError(M,info);
446 info=DSDPVecRestoreArray(b,&bb); DSDPCHKERR(info);
447 info=DSDPVecRestoreArray(x,&xx); DSDPCHKERR(info);
448 } else {
449 DSDPNoOperationError(M);
450 }
451 info=DSDPVecSetR(x,0.0);DSDPCHKERR(info);
452 info=DSDPVecSetC(x,0.0);DSDPCHKERR(info);
453 info=DSDPEventLogEnd(hsolveevent);
454 DSDPFunctionReturn(0);
455}
456
457#undef __FUNCT__
458#define __FUNCT__ "DSDPSchurMatSolve"
467 int info;
468 DSDPFunctionBegin;
469 info=DSDPSchurMatSolveM(M,b,x);DSDPCHKERR(info);
470 info=DSDPApplySMW(M,b,x);DSDPCHKERR(info);
471 info=DSDPZeroFixedVariables(M,x);DSDPCHKERR(info);
472 DSDPFunctionReturn(0);
473}
474
475#undef __FUNCT__
476#define __FUNCT__ "DSDPApplySMW"
477static int DSDPApplySMW(DSDPSchurMat M, DSDPVec rhs, DSDPVec dy){
478 int info;
479 double r=M.schur->r,rr,dr,rhsr,rssr;
480 double rhsnorm,rhsnorm3,rhs1mrhs3=0,rhs3mrhs3=0;
481 DSDPVec rhs3=M.schur->rhs3,dy3=M.schur->dy3;
482 DSDPFunctionBegin;
483
484 info=DSDPVecNormInfinity(rhs,&rhsnorm);DSDPCHKERR(info);
485 info=DSDPVecNormInfinity(rhs3,&rhsnorm3);DSDPCHKERR(info);
486 if (r==0 || rhsnorm==0){
487 info=DSDPVecSetR(dy,0); DSDPCHKERR(info);
488 info=DSDPVecSetR(rhs,0); DSDPCHKERR(info);
489 } else if (0 && rhsnorm3==0){ /* dsdp->UsePenalty==DSDPNever */
490 info=DSDPVecGetR(rhs,&rr); DSDPCHKERR(info);
491 info=DSDPVecSetR(dy,rr); DSDPCHKERR(info);
492 } else {
493 /* Use bigM penalty method and Sherman-Morrison-Woodbury */
494 info=DSDPVecGetR(rhs,&rhsr); DSDPCHKERR(info);
495 info=DSDPVecGetR(rhs3,&rssr); DSDPCHKERR(info);
496 info=DSDPVecDot(rhs3,dy,&rhs1mrhs3); DSDPCHKERR(info);
497 info=DSDPVecDot(rhs3,dy3,&rhs3mrhs3); DSDPCHKERR(info);
498 if (rssr-rhs3mrhs3==0) rssr*=(1.00001);
499 dr=-(rhs1mrhs3-rhsr )/(rssr-rhs3mrhs3);
500 info=DSDPVecAXPY(-dr,dy3,dy);DSDPCHKERR(info);
501 info=DSDPVecSetR(dy,dr); DSDPCHKERR(info);
502 info=DSDPVecSetR(rhs,rhsr); DSDPCHKERR(info);
503 info=DSDPVecDot(rhs,dy,&rhs3mrhs3); DSDPCHKERR(info);
504 if (rhs3mrhs3 <=0){
505 DSDPLogInfo(0,3,"DSDP Step Direction Not Descent, Adjusting. \n");
506 info=DSDPVecAddR(rhs3,rssr*0.1);DSDPCHKERR(info);
507 info=DSDPVecAXPY(dr,dy3,dy);DSDPCHKERR(info);
508 info=DSDPVecSetR(dy,0); DSDPCHKERR(info);
509 info=DSDPApplySMW(M,rhs,dy);DSDPCHKERR(info);
510 }
511 }
512 DSDPFunctionReturn(0);
513}
514
515#undef __FUNCT__
516#define __FUNCT__ "DSDPZeroFixedVariables"
517int DSDPZeroFixedVariables( DSDPSchurMat M, DSDPVec dy){
518 int i,info;
519 FixedVariables *fv=&M.schur->fv;
520 DSDPFunctionBegin;
521 for (i=0;i<fv->nvars;i++){
522 info=DSDPVecSetElement(dy,fv->var[i],0.0);DSDPCHKERR(info);
523 }
524 DSDPFunctionReturn(0);
525}
526
527#undef __FUNCT__
528#define __FUNCT__ "DSDPApplyFixedVariables"
529int DSDPApplyFixedVariables( DSDPSchurMat M, DSDPVec y){
530 int i,jj,info;
531 double vv,scl;
532 FixedVariables *fv=&M.schur->fv;
533 info=DSDPVecGetC(y,&scl);DSDPCHKERR(info);
534 DSDPFunctionBegin;
535 for (i=0;i<fv->nvars;i++){
536 vv=fv->fval[i]*fabs(scl);
537 jj=fv->var[i];
538 info=DSDPVecSetElement(y,jj,vv);DSDPCHKERR(info);
539 }
540 DSDPFunctionReturn(0);
541}
542
543#undef __FUNCT__
544#define __FUNCT__ "DSDPFixedVariableNorm"
545int DSDPFixedVariablesNorm( DSDPSchurMat M, DSDPVec y){
546 int i,jj,info;
547 double vv;
548 FixedVariables *fv=&M.schur->fv;
549 DSDPFunctionBegin;
550 for (i=0;i<fv->nvars;i++){
551 jj=fv->var[i]; vv=fv->fval[i];
552 info=DSDPVecAddC(y,1.0);DSDPCHKERR(info);
553 info=DSDPVecAddElement(y,jj,vv*vv);DSDPCHKERR(info);
554 }
555 DSDPFunctionReturn(0);
556}
557
558#undef __FUNCT__
559#define __FUNCT__ "DSDPComputeFixedYX"
560int DSDPComputeFixedYX( DSDPSchurMat M, DSDPVec berr){
561 int i,jj,info;
562 double vv;
563 FixedVariables *fv=&M.schur->fv;
564 DSDPFunctionBegin;
565 for (i=0;i<fv->nvars;i++){
566 jj=fv->var[i];
567 info=DSDPVecGetElement(berr,jj,&vv);DSDPCHKERR(info);
568 info=DSDPVecSetElement(berr,jj,0);DSDPCHKERR(info);
569 info=DSDPVecAddC(berr,-vv*fv->fval[i]);DSDPCHKERR(info);
570 info=DSDPVecAddR(berr,fabs(vv));DSDPCHKERR(info);
571 fv->fdual[i]=-vv;
572 if (fv->xout) fv->xout[i]=-vv;
573 DSDPLogInfo(0,2,"FIXED VAR DUAL: %d %4.4f, ADD %4.4f to objective.\n",jj,vv,-vv*fv->fval[i]);
574 }
575 DSDPFunctionReturn(0);
576}
577
578#undef __FUNCT__
579#define __FUNCT__ "DSDPIsFixed"
580int DSDPIsFixed( DSDPSchurMat M, int vari, DSDPTruth *flag){
581 int i;
582 FixedVariables *fv=&M.schur->fv;
583 DSDPFunctionBegin;
584 *flag=DSDP_FALSE;
585 for (i=0;i<fv->nvars;i++){
586 if (fv->var[i]==vari){
587 *flag=DSDP_TRUE;
588 break;
589 }
590 }
591 DSDPFunctionReturn(0);
592}
593#undef __FUNCT__
594#define __FUNCT__ "DSDPInitializeFixedVariables"
595int DSDPInitializeFixedVariable( FixedVariables *fv){
596 DSDPFunctionBegin;
597 fv->nmaxvars=0;
598 fv->nvars=0;
599 fv->fval=0;
600 fv->var=0;
601 fv->fdual=0;
602 DSDPFunctionReturn(0);
603}
604
605#undef __FUNCT__
606#define __FUNCT__ "DSDPAddFixedVariables"
607int DSDPAddFixedVariable( DSDPSchurMat M, int vari, double val){
608 int i,t,*iinew,info,nvars;
609 double *ddnew,*vvnew;
610 FixedVariables *fv=&M.schur->fv;
611 DSDPFunctionBegin;
612 nvars=fv->nvars;
613 if (nvars>=fv->nmaxvars){
614 t=2*nvars + 2;
615 DSDPCALLOC2(&iinew,int,t,&info);
616 DSDPCALLOC2(&ddnew,double,t,&info);
617 DSDPCALLOC2(&vvnew,double,t,&info);
618 for (i=0;i<nvars;i++){
619 iinew[i]=fv->var[i];
620 ddnew[i]=fv->fval[i];
621 vvnew[i]=fv->fdual[i];
622 }
623 DSDPFREE(&fv->var,&info);DSDPCHKERR(info);
624 DSDPFREE(&fv->fval,&info);DSDPCHKERR(info);
625 DSDPFREE(&fv->fdual,&info);DSDPCHKERR(info);
626 fv->var=iinew;
627 fv->fval=ddnew;
628 fv->fdual=vvnew;
629 fv->nmaxvars=t;
630 }
631 fv->var[fv->nvars]=vari;
632 fv->fval[fv->nvars]=val;
633 fv->nvars++;
634 DSDPFunctionReturn(0);
635}
636
637#include "dsdp.h"
638
639#undef __FUNCT__
640#define __FUNCT__ "DSDPSparsityInSchurMat"
649int DSDPSparsityInSchurMat(DSDP dsdp, int row, int rnnz[], int mm){
650 int info,*iptr,m=mm+2;
651 double *dd;
652 DSDPVec R=dsdp->M.schur->rhs3;
653 DSDPFunctionBegin;
654 info=DSDPVecZero(R);DSDPCHKERR(info);
655 info=DSDPVecGetArray(R,&dd);DSDPCHKERR(info);
656 iptr=(int*)dd;
657 info=DSDPSchurSparsity(dsdp,row+1,iptr,m);DSDPCHKERR(info);
658 memcpy((void*)rnnz,(void*)(iptr+1),(mm)*sizeof(int));
659 info=DSDPVecRestoreArray(R,&dd);DSDPCHKERR(info);
660 DSDPFunctionReturn(0);
661}
662
663#include "dsdp5.h"
664#undef __FUNCT__
665#define __FUNCT__ "DSDPSetFixedVariable"
675int DSDPSetFixedVariable(DSDP dsdp, int vari, double val){
676 int info;
677 DSDPFunctionBegin;
678 DSDPLogInfo(0,2,"Set Fixed Variable: %d, %12.8f\n",vari,val);
679 info= DSDPAddFixedVariable(dsdp->M,vari,val);DSDPCHKERR(info);
680 DSDPFunctionReturn(0);
681 }
682
683#undef __FUNCT__
684#define __FUNCT__ "DSDPSetFixedVariables"
695int DSDPSetFixedVariables(DSDP dsdp, double vars[], double vals[], double xout[], int nvars){
696 int i,info;
697 DSDPFunctionBegin;
698 for (i=0;i<nvars;i++){
699 info=DSDPSetFixedVariable(dsdp,(int)vars[i],vals[i]);
700 dsdp->M.schur->fv.xout=xout;
701 }
702 DSDPFunctionReturn(0);
703 }
704
705#undef __FUNCT__
706#define __FUNCT__ "DSDPGetFixedYX"
707int DSDPGetFixedYX( DSDP dsdp, int vari, double *dd){
708 int i;
709 FixedVariables *fv=&dsdp->M.schur->fv;
710 DSDPFunctionBegin;
711 for (i=0;i<fv->nvars;i++){
712 if (vari==fv->var[i]){
713 *dd=fv->fdual[i]; break;
714 }
715 }
716 DSDPFunctionReturn(0);
717}
The API to DSDP for those applications using DSDP as a subroutine library.
Internal data structure for the DSDP solver.
int DSDPSchurSparsity(DSDP, int, int[], int)
Each cone should print its state.
Definition dsdpcops.c:474
Solver, solution types, termination codes,.
DSDPTruth
Boolean variables.
@ DSDP_FALSE
@ DSDP_TRUE
int DSDPSchurMatReducePVec(DSDPSchurMat M, DSDPVec x)
Collect elements of the vector.
int DSDPSchurMatSetup(DSDPSchurMat M, DSDPVec Y)
Set up the data structure.
int DSDPSchurMatMultiply(DSDPSchurMat M, DSDPVec x, DSDPVec y)
Multiply M by a vector. y = M x.
int DSDPSchurMatSolve(DSDPSchurMat M, DSDPVec b, DSDPVec x)
Solve the linear system.
int DSDPSchurMatSetR(DSDPSchurMat M, double rr)
Set up the data structure.
int DSDPSparsityInSchurMat(DSDP dsdp, int row, int rnnz[], int mm)
Identify nonzero elements in a row of the Schur complement.
int DSDPSchurMatZeroEntries(DSDPSchurMat M)
Zero all element in the matrix.
int DSDPSchurMatInParallel(DSDPSchurMat M, DSDPTruth *flag)
Determine whether M is computed in parallel.
int DSDPSchurMatInitialize(DSDPSchurMat *M)
Initialize pointers to null.
int DSDPSchurMatDestroy(DSDPSchurMat *M)
Free the memory in the data structure.
int DSDPSchurMatAssemble(DSDPSchurMat M)
Final assembly of M.
int DSDPSchurMatFactor(DSDPSchurMat M, DSDPTruth *successful)
Factor M.
int DSDPSchurMatSetData(DSDPSchurMat *M, struct DSDPSchurMat_Ops *ops, void *data)
Set the Schur matrix with an opaque pointer and structure of function pointers.
int DSDPSchurMatView(DSDPSchurMat M)
Print the matrix.
int DSDPSchurMatOpsInitialize(struct DSDPSchurMat_Ops *dops)
Initialize function pointers to 0.
int DSDPSchurMatRowScaling(DSDPSchurMat M, DSDPVec D)
Identify which rows on on this processor.
int DSDPSchurMatShiftDiagonal(DSDPSchurMat M, double dd)
Add a scalar to each diagonal element of the matrix.
Methods of a Schur Matrix.
int DSDPSchurMatDiagonalScaling(DSDPSchurMat, DSDPVec)
Get the scaling and nonzero pattern of each diagonal element of the matrix.
Function pointers that a Schur complement matrix (dense, sparse, parallel dense) must provide.
Error handling, printing, and profiling.
struct DSDPVec_C DSDPVec
This object hold m+2 variables: a scaling of C, the y variables, and r.
Definition dsdpvec.h:25
Schur complement matrix whose solution is the Newton direction.
Internal structures for the DSDP solver.
Definition dsdp.h:65