1 /* 2 * ECOS - Embedded Conic Solver. 3 * Copyright (C) 2012-2015 A. Domahidi [domahidi@embotech.com], 4 * Automatic Control Lab, ETH Zurich & embotech GmbH, Zurich, Switzerland. 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 3 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, see <http://www.gnu.org/licenses/>. 18 */ 19 20 21 #ifndef __ECOS_H__ 22 #define __ECOS_H__ 23 24 #ifdef __cplusplus 25 extern "C" { 26 #endif 27 28 #include "glblopts.h" 29 #include "spla.h" 30 #include "cone.h" 31 #include "kkt.h" 32 33 #if PROFILING > 0 34 #include "timer.h" 35 #endif 36 37 #if CTRLC > 0 38 #include "ctrlc.h" 39 #endif 40 41 /* ECOS VERSION NUMBER - FORMAT: X.Y.Z --------------------------------- */ 42 #define ECOS_VERSION ("2.0.8") 43 44 /* DEFAULT SOLVER PARAMETERS AND SETTINGS STRUCT ----------------------- */ 45 #define MAXIT (100) /* maximum number of iterations */ 46 #define FEASTOL (1E-8) /* primal/dual infeasibility tolerance */ 47 #define ABSTOL (1E-8) /* absolute tolerance on duality gap */ 48 #define RELTOL (1E-8) /* relative tolerance on duality gap */ 49 #define FTOL_INACC (1E-4) /* inaccurate solution feasibility tol. */ 50 #define ATOL_INACC (5E-5) /* inaccurate solution absolute tol. */ 51 #define RTOL_INACC (5E-5) /* inaccurate solution relative tol. */ 52 #define GAMMA (0.99) /* scaling the final step length */ 53 #define STATICREG (1) /* static regularization: 0:off, 1:on */ 54 #define DELTASTAT (7E-8) /* regularization parameter */ 55 #define DELTA (2E-7) /* dyn. regularization parameter */ 56 #define EPS (1E-13) /* dyn. regularization threshold (do not 0!) */ 57 #define VERBOSE (1) /* bool for verbosity; PRINTLEVEL < 3 */ 58 #define NITREF (9) /* number of iterative refinement steps */ 59 #define IRERRFACT (6) /* factor by which IR should reduce err */ 60 #define LINSYSACC (1E-14) /* rel. accuracy of search direction */ 61 #define SIGMAMIN (1E-4) /* always do some centering */ 62 #define SIGMAMAX (1.0) /* never fully center */ 63 #define STEPMIN (1E-6) /* smallest step that we do take */ 64 #define STEPMAX (0.999) /* largest step allowed, also in affine dir. */ 65 #define SAFEGUARD (500) /* Maximum increase in PRES before 66 ECOS_NUMERICS is thrown. */ 67 /*Ecos exponential cone default settings*/ 68 #ifdef EXPCONE 69 #define MAX_BK (90) /*Maximum backtracking steps*/ 70 #define BK_SCALE (0.8) /*Backtracking constant*/ 71 #define MIN_DISTANCE (0.1) /* dont let sqrt(r), sqrt(-u) or sqrt(v) 72 become smaller than 73 MIN_DISTANCE*mu*/ 74 #define CENTRALITY (1) /*Centrality requirement*/ 75 #endif 76 77 78 79 /* EQUILIBRATION METHOD ------------------------------------------------ */ 80 #define EQUILIBRATE (1) /* use equlibration of data matrices? >0: yes */ 81 #define EQUIL_ITERS (3) /* number of equilibration iterations */ 82 #define RUIZ_EQUIL /* define algorithm to use - if both are ... */ 83 /*#define ALTERNATING_EQUIL*/ /* ... commented out no equlibration is used */ 84 85 86 /* EXITCODES ----------------------------------------------------------- */ 87 #define ECOS_OPTIMAL (0) /* Problem solved to optimality */ 88 #define ECOS_PINF (1) /* Found certificate of primal infeasibility */ 89 #define ECOS_DINF (2) /* Found certificate of dual infeasibility */ 90 #define ECOS_INACC_OFFSET (10) /* Offset exitflag at inaccurate results */ 91 #define ECOS_MAXIT (-1) /* Maximum number of iterations reached */ 92 #define ECOS_NUMERICS (-2) /* Search direction unreliable */ 93 #define ECOS_OUTCONE (-3) /* s or z got outside the cone, numerics? */ 94 #define ECOS_SIGINT (-4) /* solver interrupted by a signal/ctrl-c */ 95 #define ECOS_FATAL (-7) /* Unknown problem in solver */ 96 97 98 /* SETTINGS STRUCT ----------------------------------------------------- */ 99 typedef struct settings{ 100 pfloat gamma; /* scaling the final step length */ 101 pfloat delta; /* regularization parameter */ 102 pfloat eps; /* regularization threshold */ 103 pfloat feastol; /* primal/dual infeasibility tolerance */ 104 pfloat abstol; /* absolute tolerance on duality gap */ 105 pfloat reltol; /* relative tolerance on duality gap */ 106 pfloat feastol_inacc; /* primal/dual infeasibility relaxed tolerance */ 107 pfloat abstol_inacc; /* absolute relaxed tolerance on duality gap */ 108 pfloat reltol_inacc; /* relative relaxed tolerance on duality gap */ 109 idxint nitref; /* number of iterative refinement steps */ 110 idxint maxit; /* maximum number of iterations */ 111 idxint verbose; /* verbosity bool for PRINTLEVEL < 3 */ 112 #ifdef EXPCONE /*Exponential cone settings*/ 113 idxint max_bk_iter; /* Maximum backtracking iterations */ 114 pfloat bk_scale; /* Backtracking scaling */ 115 pfloat centrality; /* Centrality bound, ignored when centrality vars = 0*/ 116 #endif 117 } settings; 118 119 120 /* INFO STRUCT --------------------------------------------------------- */ 121 typedef struct stats{ 122 pfloat pcost; 123 pfloat dcost; 124 pfloat pres; 125 pfloat dres; 126 pfloat pinf; 127 pfloat dinf; 128 pfloat pinfres; 129 pfloat dinfres; 130 pfloat gap; 131 pfloat relgap; 132 pfloat sigma; 133 pfloat mu; 134 pfloat step; 135 pfloat step_aff; 136 pfloat kapovert; 137 idxint iter; 138 idxint nitref1; 139 idxint nitref2; 140 idxint nitref3; 141 #if PROFILING > 0 142 pfloat tsetup; 143 pfloat tsolve; 144 #endif 145 #if PROFILING > 1 146 pfloat tfactor; 147 pfloat tkktsolve; 148 pfloat torder; 149 pfloat tkktcreate; 150 pfloat ttranspose; 151 pfloat tperm; 152 pfloat tfactor_t1; 153 pfloat tfactor_t2; 154 #endif 155 #ifdef EXPCONE 156 /* Counters for backtracking, each of these counts 157 * one condition that can fail and cause a backtrack 158 */ 159 idxint pob; /* Potential decreases */ 160 idxint cb; /* Centrality violations */ 161 idxint cob; /* The s'z of one cone is too small w.r.t. mu */ 162 idxint pb; /* Primal infeasibility */ 163 idxint db; /* Dual infeasibility */ 164 idxint affBack; /* Total affine backtracking steps */ 165 idxint cmbBack; /* Total combined backtracking steps */ 166 167 pfloat centrality; /*Centrality at the end of the backtracking*/ 168 #endif 169 170 } stats; 171 172 173 /* ALL DATA NEEDED BY SOLVER ------------------------------------------- */ 174 typedef struct pwork{ 175 /* dimensions */ 176 idxint n; /* number of primal variables x */ 177 idxint m; /* number of conically constrained variables s */ 178 idxint p; /* number of equality constraints */ 179 idxint D; /* degree of the cone */ 180 181 /* variables */ 182 pfloat* x; /* primal variables */ 183 pfloat* y; /* multipliers for equality constaints */ 184 pfloat* z; /* multipliers for conic inequalities */ 185 pfloat* s; /* slacks for conic inequalities */ 186 pfloat* lambda; /* scaled variable */ 187 pfloat kap; /* kappa (homogeneous embedding) */ 188 pfloat tau; /* tau (homogeneous embedding) */ 189 190 /* best iterate seen so far */ 191 /* variables */ 192 pfloat* best_x; /* primal variables */ 193 pfloat* best_y; /* multipliers for equality constaints */ 194 pfloat* best_z; /* multipliers for conic inequalities */ 195 pfloat* best_s; /* slacks for conic inequalities */ 196 pfloat best_kap; /* kappa (homogeneous embedding) */ 197 pfloat best_tau; /* tau (homogeneous embedding) */ 198 pfloat best_cx; 199 pfloat best_by; 200 pfloat best_hz; 201 stats* best_info; /* info of best iterate */ 202 203 /* temporary stuff holding search direction etc. */ 204 pfloat* dsaff; 205 pfloat* dzaff; 206 pfloat* W_times_dzaff; 207 pfloat* dsaff_by_W; 208 pfloat* saff; 209 pfloat* zaff; 210 211 /* cone */ 212 cone* C; 213 214 /* problem data */ 215 spmat* A; spmat* G; pfloat* c; pfloat* b; pfloat* h; 216 217 /* indices that map entries of A and G to the KKT matrix */ 218 idxint *AtoK; idxint *GtoK; 219 220 #if defined EQUILIBRATE && EQUILIBRATE > 0 221 /* equilibration vector */ 222 pfloat *xequil; 223 pfloat *Aequil; 224 pfloat *Gequil; 225 #endif 226 227 /* scalings of problem data */ 228 pfloat resx0; pfloat resy0; pfloat resz0; 229 230 /* residuals */ 231 pfloat *rx; pfloat *ry; pfloat *rz; pfloat rt; 232 pfloat hresx; pfloat hresy; pfloat hresz; 233 234 /* norm iterates */ 235 pfloat nx,ny,nz,ns; 236 237 /* temporary storage */ 238 pfloat cx; pfloat by; pfloat hz; pfloat sz; 239 240 /* KKT System */ 241 kkt* KKT; 242 243 /* info struct */ 244 stats* info; 245 246 /* settings struct */ 247 settings* stgs; 248 249 } pwork; 250 251 252 /* SOME USEFUL MACROS -------------------------------------------------- */ 253 #define MAX(X,Y) ((X) < (Y) ? (Y) : (X)) /* maximum of 2 expressions */ 254 /* safe division x/y where y is assumed to be positive! */ 255 #define SAFEDIV_POS(X,Y) ( (Y) < EPS ? ((X)/EPS) : (X)/(Y) ) 256 257 258 /* METHODS */ 259 260 /* set up work space 261 * could be done by codegen 262 * 263 * Parameters: 264 * idxint n Number of variables 265 * idxint m Number of inequalities, number of rows of G 266 * idxint p Number of equality constraints 267 * idxint l Dimension of positive orthant 268 * idxint ncones Number of second order cones 269 * idxint* q Array of length 'ncones', defines the dimension of each cone 270 * idxint nex Number of exponential cones 271 * pfloat* Gpr Sparse G matrix data array (column compressed storage) 272 * idxint* Gjc Sparse G matrix column index array (column compressed storage) 273 * idxint* Gir Sparse G matrix row index array (column compressed storage) 274 * pfloat* Apr Sparse A matrix data array (column compressed storage) (can be all NULL if no equalities are present) 275 * idxint* Ajc Sparse A matrix column index array (column compressed storage) (can be all NULL if no equalities are present) 276 * idxint* Air Sparse A matrix row index array (column compressed storage) (can be all NULL if no equalities are present) 277 * pfloat* c Array of size n, cost function weights 278 * pfloat* h Array of size m, RHS vector of cone constraint 279 * pfloat* b Array of size p, RHS vector of equalities (can be NULL if no equalities are present) 280 */ 281 pwork* ECOS_setup(idxint n, idxint m, idxint p, idxint l, idxint ncones, idxint* q, idxint nex, 282 pfloat* Gpr, idxint* Gjc, idxint* Gir, 283 pfloat* Apr, idxint* Ajc, idxint* Air, 284 pfloat* c, pfloat* h, pfloat* b); 285 286 287 #ifdef EXPCONE 288 pfloat expConeLineSearch(pwork* w, pfloat dtau, pfloat dkappa, idxint affine); 289 #endif 290 291 /* solve */ 292 idxint ECOS_solve(pwork* w); 293 294 /** 295 * Cleanup: free memory (not used for embedded solvers, only standalone) 296 * 297 * Use the second argument to give the number of variables to NOT free. 298 * This is useful if you want to use the result of the optimization without 299 * copying over the arrays. One use case is the MEX interface, where we 300 * do not want to free x,y,s,z (depending on the number of LHS). 301 */ 302 void ECOS_cleanup(pwork* w, idxint keepvars); 303 304 305 /** 306 * Version: returns the current version number 307 * Use a character array of length 7 to obtain the version number 308 * in the format 309 * x.y.zzz 310 * where x is the major, y the minor and zzz the build number 311 */ 312 const char* ECOS_ver(void); 313 314 315 /* ------------------- EXPERT LEVEL INTERFACES ---------------------- */ 316 317 /* 318 * Updates one element of the RHS vector h of inequalities 319 * After the call, w->h[idx] = value (but equilibrated) 320 */ 321 void ecos_updateDataEntry_h(pwork* w, idxint idx, pfloat value); 322 323 /* 324 * Updates one element of the OBJ vector c of inequalities 325 * After the call, w->c[idx] = value (but equilibrated) 326 */ 327 void ecos_updateDataEntry_c(pwork* w, idxint idx, pfloat value); 328 329 /* 330 * Updates numerical data for G, A, c, h, and b, 331 * and re-equilibrates. 332 * Then updates the corresponding KKT entries. 333 */ 334 void ECOS_updateData(pwork *w, pfloat *Gpr, pfloat *Apr, 335 pfloat* c, pfloat* h, pfloat* b); 336 337 338 #ifdef __cplusplus 339 } 340 #endif 341 342 #endif 343 344