11 static int hfactorevent=0,hsolveevent=0;
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);} } 20 #define __FUNCT__ "DSDPSchurMatSetData" 28 int DSDPSchurMatSetData(
DSDPSchurMat *M,
struct DSDPSchurMat_Ops* ops,
void*data){
32 DSDPFunctionReturn(0);
35 static const char* schurmatname=
"NOT NAMED YET";
38 #define __FUNCT__ "DSDPSchurMatOpsInitialize" 46 if (dops==NULL)
return 0;
48 dops->matrownonzeros=0;
50 dops->mataddelement=0;
51 dops->matadddiagonal=0;
52 dops->matshiftdiagonal=0;
54 dops->matscaledmultiply=0;
58 dops->pmatonprocessor=0;
59 dops->pmatwhichdiag=0;
60 dops->pmatdistributed=0;
65 dops->matname=schurmatname;
66 DSDPFunctionReturn(0);
69 static struct DSDPSchurMat_Ops dsdpmops;
73 #define __FUNCT__ "DSDPSchurMatOpsInitialize" 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);
91 #define __FUNCT__ "DSDPSchurMatZeroEntries" 100 if (M.dsdpops->matzero){
101 info=(M.dsdpops->matzero)(M.data); DSDPChkMatError(M,info);
103 DSDPNoOperationError(M);
105 DSDPFunctionReturn(0);
110 #define __FUNCT__ "DSDPSchurMatShiftDiagonal" 123 if (dd==0){DSDPFunctionReturn(0);}
125 if (M.dsdpops->matshiftdiagonal){
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);
130 DSDPNoOperationError(M);
132 DSDPFunctionReturn(0);
137 #define __FUNCT__ "DSDPSchurMatInParallel" 152 if (M.dsdpops->pmatdistributed){
153 info=(M.dsdpops->pmatdistributed)(M.data,&flg); DSDPChkMatError(M,info);
158 DSDPFunctionReturn(0);
164 #define __FUNCT__ "DSDPSchurMatAssemble" 177 if (M.dsdpops->matassemble){
178 info=(M.dsdpops->matassemble)(M.data); DSDPChkMatError(M,info);
180 DSDPNoOperationError(M);
182 DSDPFunctionReturn(0);
186 #define __FUNCT__ "DSDPSchurMatFactor" 198 DSDPVec rhs3=M.schur->rhs3,dy3=M.schur->dy3;
201 DSDPEventLogBegin(hfactorevent);
202 if (M.dsdpops->matfactor){
203 info=(M.dsdpops->matfactor)(M.data,&flag); DSDPChkMatError(M,info);
206 DSDPLogInfo(0,2,
"Indefinite Schur Matrix -- Bad Factorization\n");
209 DSDPNoOperationError(M);
211 DSDPEventLogEnd(hfactorevent);
213 info=DSDPSchurMatSolveM(M,rhs3,dy3);DSDPCHKERR(info);}
214 else {info=DSDPVecZero(dy3);DSDPCHKERR(info);}
215 DSDPFunctionReturn(0);
220 #define __FUNCT__ "DSDPSchurMatMultiply" 233 double *xx,*yy,r=M.schur->r;
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);
245 info=DSDPVecRestoreArray(y,&yy); DSDPCHKERR(info);
246 info=DSDPVecRestoreArray(x,&xx); DSDPCHKERR(info);
248 DSDPNoOperationError(M);
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);
258 DSDPFunctionReturn(0);
262 #define __FUNCT__ "DSDPSchurMatMultR" 265 double *xx,*yy,r=M.schur->r;
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);
277 info=DSDPVecRestoreArray(y,&yy); DSDPCHKERR(info);
278 info=DSDPVecRestoreArray(x,&xx); DSDPCHKERR(info);
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);
288 info=DSDPVecZero(y);DSDPCHKERR(info);
291 DSDPFunctionReturn(0);
295 #define __FUNCT__ "DSDPSchurMatReducePVec" 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);
321 DSDPNoOperationError(M);
324 info=DSDPZeroFixedVariables(M,x);DSDPCHKERR(info);
325 DSDPFunctionReturn(0);
331 #define __FUNCT__ "DSDPSchurMatSetR" 341 DSDPFunctionReturn(0);
345 #define __FUNCT__ "DSDPSchurMatSetup" 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);
361 DSDPNoOperationError(M);
363 DSDPEventLogRegister(
"Factor Newton Eq.",&hfactorevent);
364 DSDPEventLogRegister(
"Solve Newton Eq.",&hsolveevent);
365 DSDPFunctionReturn(0);
370 #define __FUNCT__ "DSDPSchurMatView" 379 if (M.dsdpops->matview){
380 info=(M.dsdpops->matview)(M.data); DSDPChkMatError(M,info);
382 DSDPNoOperationError(M);
384 info=DSDPVecView(M.schur->rhs3);DSDPCHKERR(info);
385 DSDPFunctionReturn(0);
390 #define __FUNCT__ "DSDPSchurMatRowScaling" 403 info=DSDPZeroFixedVariables(M,D);DSDPCHKERR(info);
404 DSDPFunctionReturn(0);
408 #define __FUNCT__ "DSDPSchurMatDestroy" 417 if ((*M).dsdpops->matdestroy){
418 info=((*M).dsdpops->matdestroy)((*M).data); DSDPChkMatError(*M,info);
424 info=DSDPVecDestroy(&M->schur->rhs3);DSDPCHKERR(info);
425 info=DSDPVecDestroy(&M->schur->dy3);DSDPCHKERR(info);
428 info=DSDPSchurMatSetData(M,&dsdpmops,0); DSDPCHKERR(info);
429 DSDPFREE(&M->schur,&info);DSDPCHKERR(info);
430 DSDPFunctionReturn(0);
434 #define __FUNCT__ "DSDPSchurMatSolveM" 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);
449 DSDPNoOperationError(M);
451 info=DSDPVecSetR(x,0.0);DSDPCHKERR(info);
452 info=DSDPVecSetC(x,0.0);DSDPCHKERR(info);
453 info=DSDPEventLogEnd(hsolveevent);
454 DSDPFunctionReturn(0);
458 #define __FUNCT__ "DSDPSchurMatSolve" 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);
476 #define __FUNCT__ "DSDPApplySMW" 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;
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){
490 info=DSDPVecGetR(rhs,&rr); DSDPCHKERR(info);
491 info=DSDPVecSetR(dy,rr); DSDPCHKERR(info);
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);
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);
512 DSDPFunctionReturn(0);
516 #define __FUNCT__ "DSDPZeroFixedVariables" 519 FixedVariables *fv=&M.schur->fv;
521 for (i=0;i<fv->nvars;i++){
522 info=DSDPVecSetElement(dy,fv->var[i],0.0);DSDPCHKERR(info);
524 DSDPFunctionReturn(0);
528 #define __FUNCT__ "DSDPApplyFixedVariables" 532 FixedVariables *fv=&M.schur->fv;
533 info=DSDPVecGetC(y,&scl);DSDPCHKERR(info);
535 for (i=0;i<fv->nvars;i++){
536 vv=fv->fval[i]*fabs(scl);
538 info=DSDPVecSetElement(y,jj,vv);DSDPCHKERR(info);
540 DSDPFunctionReturn(0);
544 #define __FUNCT__ "DSDPFixedVariableNorm" 548 FixedVariables *fv=&M.schur->fv;
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);
555 DSDPFunctionReturn(0);
559 #define __FUNCT__ "DSDPComputeFixedYX" 563 FixedVariables *fv=&M.schur->fv;
565 for (i=0;i<fv->nvars;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);
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]);
575 DSDPFunctionReturn(0);
579 #define __FUNCT__ "DSDPIsFixed" 582 FixedVariables *fv=&M.schur->fv;
585 for (i=0;i<fv->nvars;i++){
586 if (fv->var[i]==vari){
591 DSDPFunctionReturn(0);
594 #define __FUNCT__ "DSDPInitializeFixedVariables" 595 int DSDPInitializeFixedVariable( FixedVariables *fv){
602 DSDPFunctionReturn(0);
606 #define __FUNCT__ "DSDPAddFixedVariables" 607 int DSDPAddFixedVariable(
DSDPSchurMat M,
int vari,
double val){
608 int i,t,*iinew,info,nvars;
609 double *ddnew,*vvnew;
610 FixedVariables *fv=&M.schur->fv;
613 if (nvars>=fv->nmaxvars){
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++){
620 ddnew[i]=fv->fval[i];
621 vvnew[i]=fv->fdual[i];
623 DSDPFREE(&fv->var,&info);DSDPCHKERR(info);
624 DSDPFREE(&fv->fval,&info);DSDPCHKERR(info);
625 DSDPFREE(&fv->fdual,&info);DSDPCHKERR(info);
631 fv->var[fv->nvars]=vari;
632 fv->fval[fv->nvars]=val;
634 DSDPFunctionReturn(0);
640 #define __FUNCT__ "DSDPSparsityInSchurMat" 650 int info,*iptr,m=mm+2;
654 info=DSDPVecZero(R);DSDPCHKERR(info);
655 info=DSDPVecGetArray(R,&dd);DSDPCHKERR(info);
658 memcpy((
void*)rnnz,(
void*)(iptr+1),(mm)*
sizeof(
int));
659 info=DSDPVecRestoreArray(R,&dd);DSDPCHKERR(info);
660 DSDPFunctionReturn(0);
665 #define __FUNCT__ "DSDPSetFixedVariable" 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);
684 #define __FUNCT__ "DSDPSetFixedVariables" 698 for (i=0;i<nvars;i++){
700 dsdp->M.schur->fv.xout=xout;
702 DSDPFunctionReturn(0);
706 #define __FUNCT__ "DSDPGetFixedYX" 707 int DSDPGetFixedYX(
DSDP dsdp,
int vari,
double *dd){
709 FixedVariables *fv=&dsdp->M.schur->fv;
711 for (i=0;i<fv->nvars;i++){
712 if (vari==fv->var[i]){
713 *dd=fv->fdual[i];
break;
716 DSDPFunctionReturn(0);
DSDPTruth
Boolean variables.
struct DSDPVec_C DSDPVec
This object hold m+2 variables: a scaling of C, the y variables, and r.
Schur complement matrix whose solution is the Newton direction.
int DSDPSchurMatSolve(DSDPSchurMat M, DSDPVec b, DSDPVec x)
Solve the linear system.
int DSDPSchurMatFactor(DSDPSchurMat M, DSDPTruth *successful)
Factor M.
Error handling, printing, and profiling.
int DSDPSparsityInSchurMat(DSDP dsdp, int row, int rnnz[], int mm)
Identify nonzero elements in a row of the Schur complement.
Function pointers that a Schur complement matrix (dense, sparse, parallel dense) must provide.
Internal structures for the DSDP solver.
int DSDPSchurMatMultiply(DSDPSchurMat M, DSDPVec x, DSDPVec y)
Multiply M by a vector. y = M x.
The API to DSDP for those applications using DSDP as a subroutine library.
int DSDPSchurMatZeroEntries(DSDPSchurMat M)
Zero all element in the matrix.
Solver, solution types, termination codes,.
Internal data structure for the DSDP solver.
int DSDPSchurMatRowScaling(DSDPSchurMat M, DSDPVec D)
Identify which rows on on this processor.
int DSDPSchurMatView(DSDPSchurMat M)
Print the matrix.
Methods of a Schur Matrix.
int DSDPSetFixedVariables(DSDP dsdp, double vars[], double vals[], double xout[], int nvars)
Fix variable y to exact values.
int DSDPSchurMatDiagonalScaling(DSDPSchurMat, DSDPVec)
Get the scaling and nonzero pattern of each diagonal element of the matrix.
int DSDPSchurMatOpsInitialize(struct DSDPSchurMat_Ops *dops)
Initialize function pointers to 0.
int DSDPSchurMatAssemble(DSDPSchurMat M)
Final assembly of M.
int DSDPSchurMatInParallel(DSDPSchurMat M, DSDPTruth *flag)
Determine whether M is computed in parallel.
int DSDPSchurMatSetR(DSDPSchurMat M, double rr)
Set up the data structure.
int DSDPSchurMatSetup(DSDPSchurMat M, DSDPVec Y)
Set up the data structure.
int DSDPSchurMatReducePVec(DSDPSchurMat M, DSDPVec x)
Collect elements of the vector.
int DSDPSchurMatDestroy(DSDPSchurMat *M)
Free the memory in the data structure.
int DSDPSchurMatShiftDiagonal(DSDPSchurMat M, double dd)
Add a scalar to each diagonal element of the matrix.
int DSDPSchurSparsity(DSDP, int, int[], int)
Each cone should print its state.
int DSDPSetFixedVariable(DSDP dsdp, int vari, double val)
Fix variable y to exact value.
int DSDPSchurMatInitialize(DSDPSchurMat *M)
Initialize pointers to null.