1 /* -------------------------------------------------------------
2 
3 This file is a component of SDPA
4 Copyright (C) 2004-2012 SDPA Project
5 
6 This program is free software; you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation; either version 2 of the License, or
9 (at your option) any later version.
10 
11 This program is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
14 GNU General Public License for more details.
15 
16 You should have received a copy of the GNU General Public License
17 along with this program; if not, write to the Free Software
18 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
19 
20 ------------------------------------------------------------- */
21 // Here is the start of ``example1.cpp''
22 
23 #include <cstdio>
24 #include <cstdlib>
25 #include <sdpa_call.h>
26 
27 /*
28 example1.dat:
29 
30 "Example 1: mDim = 3, nBLOCK = 1, {2}"
31    3  =  mDIM
32    1  =  nBLOCK
33    2  = bLOCKsTRUCT
34 {48, -8, 20}
35 { {-11,  0}, { 0, 23} }
36 { { 10,  4}, { 4,  0} }
37 { {  0,  0}, { 0, -8} }
38 { {  0, -8}, {-8, -2} }
39 */
40 
41 void printVector(double* ele, int dim, char* printFormat,
42 		 FILE* fpout);
43 void printMatrix(double* ele, int dim, char* printFormat,
44 		 FILE* fpout);
45 void printDimacsError(double dimacs_error[7],char* printFormat,
46 		      FILE* fpout);
47 
main()48 int main ()
49 {
50   SDPA::printSDPAVersion(stdout);
51   SDPA	Problem1;
52   Problem1.setDisplay(stdout);
53 
54   // All parameteres are renewed
55   Problem1.setParameterType(SDPA::PARAMETER_DEFAULT);
56 
57   // If necessary, each parameter can be set independently
58 
59   // Problem1.setParameterMaxIteration(100);
60   // Problem1.setParameterEpsilonStar(1.0e-7);
61   // Problem1.setParameterLambdaStar(1.0e+2);
62   // Problem1.setParameterOmegaStar(2.0);
63   // Problem1.setParameterLowerBound(-1.0e+5);
64   // Problem1.setParameterUppwerBound(1.0e+5);
65   // Problem1.setParameterBetaStar(0.1);
66   // Problem1.setParameterBetaBar(0.2);
67   // Problem1.setParameterGammaStar(0.9);
68   // Problem1.setParameterEpsilonDash(1.0e-7);
69   // Problem1.setParameterPrintXVec((char*)"%+8.3e" );
70   // Problem1.setParameterPrintXMat((char*)"%+8.3e" );
71   // Problem1.setParameterPrintYMat((char*)"%+8.3e" );
72   // Problem1.setParameterPrintInformation((char*)"%+10.16e");
73 
74   Problem1.printParameters(stdout);
75 
76   int mDIM   = 3;
77   int nBlock = 1;
78   Problem1.inputConstraintNumber(mDIM);
79   Problem1.inputBlockNumber(nBlock);
80   Problem1.inputBlockSize(1,2);
81   Problem1.inputBlockType(1,SDPA::SDP);
82 
83   Problem1.initializeUpperTriangleSpace();
84 
85   Problem1.inputCVec(1,48);
86   Problem1.inputCVec(2,-8);
87   Problem1.inputCVec(3,20);
88 
89   Problem1.inputElement(0, 1, 1, 1, -11);
90   Problem1.inputElement(0, 1, 2, 2,  23);
91 
92   Problem1.inputElement(1, 1, 1, 1,  10);
93   Problem1.inputElement(1, 1, 1, 2,   4);
94 
95   Problem1.inputElement(2, 1, 2, 2,  -8);
96 
97   Problem1.inputElement(3, 1, 1, 2,  -8);
98   Problem1.inputElement(3, 1, 2, 2,  -2);
99 
100   Problem1.initializeUpperTriangle();
101   Problem1.initializeSolve();
102 
103   // if necessary, dump input data and initial point
104   // Problem1.writeInputSparse((char*)"tmp.dat-s",(char*)"%+8.3e");
105   // Problem1.writeInitSparse((char*)"tmp.ini-s",(char*)"%+8.3e");
106 
107   Problem1.solve();
108 
109   fprintf(stdout, "\nStop iteration = %d\n",
110 	  Problem1.getIteration());
111   char phase_string[30];
112   Problem1.getPhaseString(phase_string);
113   fprintf(stdout, "Phase          = %s\n", phase_string);
114   fprintf(stdout, "objValPrimal   = %+10.6e\n",
115 	  Problem1.getPrimalObj());
116   fprintf(stdout, "objValDual     = %+10.6e\n",
117 	  Problem1.getDualObj());
118   fprintf(stdout, "p. feas. error = %+10.6e\n",
119 	  Problem1.getPrimalError());
120   fprintf(stdout, "d. feas. error = %+10.6e\n\n",
121 	  Problem1.getDualError());
122 
123 
124   fprintf(stdout, "xVec = \n");
125   // Problem1.printResultXVec();
126   printVector(Problem1.getResultXVec(),
127 	      Problem1.getConstraintNumber(), (char*)"%+8.3e",
128 	      stdout);
129 
130   fprintf(stdout, "xMat = \n");
131   // Problem1.printResultXMat();
132   for (int l=0; l<Problem1.getBlockNumber(); ++l) {
133     if (Problem1.getBlockType(l+1) == SDPA::SDP) {
134       printMatrix(Problem1.getResultXMat(l+1),
135 		  Problem1.getBlockSize(l+1), (char*)"%+8.3e",
136 		  stdout);
137     }
138     else if (Problem1.getBlockType(l+1) == SDPA::SOCP) {
139       printf("current version does not support SOCP\n");
140     }
141     if (Problem1.getBlockType(l+1) == SDPA::LP) {
142       printVector(Problem1.getResultXMat(l+1),
143 		  Problem1.getBlockSize(l+1), (char*)"%+8.3e",
144 		  stdout);
145     }
146   }
147 
148   fprintf(stdout, "yMat = \n");
149   // Problem1.printResultYMat();
150   for (int l=0; l<Problem1.getBlockNumber(); ++l) {
151     if (Problem1.getBlockType(l+1) == SDPA::SDP) {
152       printMatrix(Problem1.getResultYMat(l+1),
153 		  Problem1.getBlockSize(l+1), (char*)"%+8.3e",
154 		  stdout);
155     }
156     else if (Problem1.getBlockType(l+1) == SDPA::SOCP) {
157       printf("current version does not support SOCP\n");
158     }
159     if (Problem1.getBlockType(l+1) == SDPA::LP) {
160       printVector(Problem1.getResultYMat(l+1),
161 		  Problem1.getBlockSize(l+1), (char*)"%+8.3e",
162 		  stdout);
163     }
164   }
165 
166   double dimacs_error[7];
167   Problem1.getDimacsError(dimacs_error);
168   printDimacsError(dimacs_error,(char*)"%+8.3e",stdout);
169 
170   // Problem1.printComputationTime(stdout);
171 
172   Problem1.terminate();
173   exit(0);
174 };
175 
printVector(double * ele,int dim,char * printFormat,FILE * fpout)176 void printVector(double* ele, int dim, char* printFormat, FILE* fpout)
177 {
178   fprintf(fpout,"[ ");
179   for (int k=0; k<dim-1; ++k) {
180     fprintf(fpout,printFormat,ele[k]);
181     fprintf(fpout," ");
182   }
183   fprintf(fpout,printFormat,ele[dim-1]);
184   fprintf(fpout,"]; \n");
185 }
186 
printMatrix(double * ele,int dim,char * printFormat,FILE * fpout)187 void printMatrix(double* ele, int dim, char* printFormat, FILE* fpout)
188 {
189   fprintf(fpout,"[\n");
190   for (int i=0; i<dim; ++i) {
191     fprintf(fpout,"[ ");
192     for (int j=0; j<dim-1; ++j) {
193       fprintf(fpout,printFormat,ele[i+dim*j]);
194       fprintf(fpout," ");
195     }
196     fprintf(fpout,printFormat,ele[i+dim*(dim-1)]);
197     fprintf(fpout,"]; \n");
198   }
199   fprintf(fpout,"]; \n");
200 }
201 
printDimacsError(double dimacs_error[7],char * printFormat,FILE * fpout)202 void printDimacsError(double dimacs_error[7],char* printFormat,
203 		      FILE* fpout)
204 {
205   fprintf(fpout,  "\n");
206   fprintf(fpout,  "* DIMACS_ERRORS * \n");
207   fprintf(fpout,  "err1 = ");
208   fprintf(fpout,  printFormat, dimacs_error[1]);
209   fprintf(fpout, "  [||Ax-b|| / (1+||b||_1)]\n");
210   fprintf(fpout,  "err2 = ");
211   fprintf(fpout,  printFormat, dimacs_error[2]);
212   fprintf(fpout, "  [max(0, -lambda(x)/(1+||b||_1))]\n");
213   fprintf(fpout,  "err3 = ");
214   fprintf(fpout,  printFormat, dimacs_error[3]);
215   fprintf(fpout, "  [||A^Ty + z - c || / (1+||c||_1)]\n");
216   fprintf(fpout,  "err4 = ");
217   fprintf(fpout,  printFormat, dimacs_error[4]);
218   fprintf(fpout, "  [max(0, -lambda(z)/(1+||c||_1))]\n");
219   fprintf(fpout,  "err5 = ");
220   fprintf(fpout,  printFormat, dimacs_error[5]);
221   fprintf(fpout, "  [(<c,x> - <b,y>) / (1 + |<c,x>| + |<b,y>|)]\n");
222   fprintf(fpout,  "err6 = ");
223   fprintf(fpout,  printFormat, dimacs_error[6]);
224   fprintf(fpout, "  [<x,z> / (1 + |<c,x>| + |<b,y>|)]\n");
225   fprintf(fpout,  "\n");
226 }
227 
228