make unique types/entities part of irprog
[libfirm] / ir / lpp / lpp_gurobi.c
1 /*
2  * Copyright (C) 2005-2011 University of Karlsruhe.  All right reserved.
3  *
4  * This file is part of libFirm.
5  *
6  * This file may be distributed and/or modified under the terms of the
7  * GNU General Public License version 2 as published by the Free Software
8  * Foundation and appearing in the file LICENSE.GPL included in the
9  * packaging of this file.
10  *
11  * Licensees holding valid libFirm Professional Edition licenses may use
12  * this file in accordance with the libFirm Commercial License.
13  * Agreement provided with the Software.
14  *
15  * This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
16  * WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR
17  * PURPOSE.
18  */
19
20 /**
21  * @file
22  * @author  Matthias Braun
23  */
24 #include "config.h"
25
26 #ifdef WITH_GUROBI
27 #include "lpp_gurobi.h"
28
29 #include <stdio.h>
30 #include <stdlib.h>
31 #include <math.h>
32
33 #include "obst.h"
34
35 #include <gurobi_c.h>
36
37 #include "error.h"
38 #include "sp_matrix.h"
39
40 static char gurobi_cst_encoding[4] = { 0, GRB_EQUAL, GRB_LESS_EQUAL, GRB_GREATER_EQUAL };
41 static char gurobi_var_encoding[4] = { 0, 0, GRB_CONTINUOUS, GRB_BINARY };
42
43 typedef struct _gurobi_t {
44         lpp_t *lpp;
45         GRBenv *env;
46         GRBenv *modelenv;
47         GRBmodel *model;
48 } gurobi_t;
49
50 static void check_gurobi_error(gurobi_t *grb, int error)
51 {
52         if (error != 0) {
53                 panic("gurobi error: %s", GRBgeterrormsg(grb->env));
54         }
55 }
56
57 static gurobi_t *new_gurobi(lpp_t *lpp)
58 {
59         int error;
60
61         gurobi_t *grb = XMALLOCZ(gurobi_t);
62         grb->lpp = lpp;
63         /* /tmp/firm_gurobi.log is a hack (see below) */
64         error = GRBloadenv(&grb->env, "/tmp/firm_gurobi.log");
65         check_gurobi_error(grb, error);
66         /* Matze: do not set the FILE* for logging output. Because:
67          *  a) the function is deprecated
68          *  b) gurobi closes the FILE handle when it is done, which leads to
69          *     very unexpected effects when you pass stdout or stderr as logging
70          *     output.
71          * The only thing gurobi sanely supports is giving a string with a filename
72          * :-( ...so we use /tmp/firm_gurobi.log as a temporary measure...
73          */
74 #if 0
75         error = GRBsetlogfile(grb->env, lpp->log);
76         check_gurobi_error(grb, error);
77 #else
78         if (lpp->log != stdout && lpp->log != stderr) {
79                 error = GRBsetintparam(grb->env, GRB_INT_PAR_OUTPUTFLAG, 0);
80                 check_gurobi_error(grb, error);
81         }
82 #endif
83
84         return grb;
85 }
86
87 static void free_gurobi(gurobi_t *grb)
88 {
89         GRBfreemodel(grb->model);
90         GRBfreeenv(grb->env);
91         free(grb);
92 }
93
94 /**
95  * Build CPLEX data structure from LPP matrix.
96  * @note: The LPP matrix is freed after this step, to save memory.
97  */
98 static void gurobi_construct(gurobi_t *grb)
99 {
100         const matrix_elem_t *elem;
101         int                  i, o;
102         //int                  sv_cnt;
103         //int                 *indices;
104         //double              *startv;
105         int                  numcols, numrows, numentries;
106         int                  objsen, *matbeg, *matcnt, *matind;
107         double               *obj, *rhs, *matval, *lb;
108         char                 *sense, *vartype;
109         char                 **colname, **rowname;
110         struct obstack       obst;
111         lpp_t                *lpp = grb->lpp;
112         int                  error;
113
114         numcols    = lpp->var_next-1;
115         numrows    = lpp->cst_next-1;
116         numentries = matrix_get_entries(lpp->m);
117         objsen     = lpp->opt_type == lpp_minimize ? 1 : -1;
118         obstack_init(&obst);
119
120         obj     = obstack_alloc(&obst, numcols * sizeof(*obj));
121         lb      = obstack_alloc(&obst, numcols * sizeof(*lb));
122         colname = obstack_alloc(&obst, numcols * sizeof(*colname));
123         rowname = obstack_alloc(&obst, numrows * sizeof(*rowname));
124         vartype = obstack_alloc(&obst, numcols * sizeof(*vartype));
125         //indices = obstack_alloc(&obst, numcols * sizeof(*indices));
126         //startv  = obstack_alloc(&obst, numcols * sizeof(*startv));
127         matbeg  = obstack_alloc(&obst, numcols * sizeof(*matbeg));
128         matcnt  = obstack_alloc(&obst, numcols * sizeof(*matcnt));
129         matind  = obstack_alloc(&obst, numentries * sizeof(*matind));
130         matval  = obstack_alloc(&obst, numentries * sizeof(*matval));
131         rhs     = obstack_alloc(&obst, numrows * sizeof(*rhs));
132         sense   = obstack_alloc(&obst, numrows * sizeof(*sense));
133
134         o      = 0;
135         //sv_cnt = 0;
136         /* fill the CPLEX matrix*/
137         for (i = 0; i < numcols; ++i) {
138                 lpp_name_t *curr_var = lpp->vars[1+i];
139
140                 obj[i] = matrix_get(lpp->m, 0, 1+i);
141                 lb[i]  = 0.0;
142
143                 colname[i] = (char*) curr_var->name;
144                 vartype[i] = gurobi_var_encoding[curr_var->type.var_type];
145
146 #if 0
147                 if (curr_var->value_kind == lpp_value_start) {
148                         panic("start values not supported in gurobi yet");
149                         indices[sv_cnt]  = i;
150                         startv[sv_cnt++] = curr_var->value;
151                 }
152 #endif
153
154                 matbeg[i] = o;
155                 matcnt[i] = 0;
156                 matrix_foreach_in_col(lpp->m, 1 + i, elem) {
157                         if (elem->row == 0)
158                                 continue;
159                         matind[o] = elem->row-1;
160                         matval[o] = elem->val;
161                         matcnt[i]++;
162                         o++;
163                 }
164         }
165
166         /* get constraint stuff (right hand side, type, name) */
167         for (i = 0; i < numrows; ++i) {
168                 lpp_name_t *curr_cst = lpp->csts[1 + i];
169
170                 rhs[i]     = matrix_get(lpp->m, 1 + i, 0);
171                 sense[i]   = gurobi_cst_encoding[curr_cst->type.cst_type];
172                 rowname[i] = (char*) curr_cst->name;
173         }
174
175         error = GRBloadmodel(grb->env, &grb->model, lpp->name, numcols, numrows,
176                              objsen, 0, obj, sense, rhs, matbeg, matcnt, matind,
177                              matval, lb, NULL, vartype, colname, rowname);
178         check_gurobi_error(grb, error);
179         grb->modelenv = GRBgetenv(grb->model);
180
181         obstack_free(&obst, NULL);
182         lpp_free_matrix(lpp);
183 }
184
185 static void gurobi_solve(gurobi_t *grb)
186 {
187         lpp_t *lpp = grb->lpp;
188         int i;
189         int optimstatus;
190         int error;
191         int numcols = lpp->var_next-1;
192         double *values;
193         double  iterations;
194
195         /* Set the time limit appropriately */
196         if(lpp->time_limit_secs > 0.0) {
197                 error = GRBsetdblparam(grb->modelenv, GRB_DBL_PAR_TIMELIMIT, lpp->time_limit_secs);
198                 check_gurobi_error(grb, error);
199         }
200
201 #if 0
202         /*
203          * If a bound of the objective function is supplied,
204          * set it accordingly, dependign on minimization or maximization.
205          */
206         if(lpp->set_bound) {
207                 fprintf(stderr, "Warning: gurobi bound not implemented yet\n");
208         }
209 #endif
210
211         /* solve */
212         error = GRBoptimize(grb->model);
213         check_gurobi_error(grb, error);
214
215         /* get solution status */
216         error = GRBgetintattr(grb->model, GRB_INT_ATTR_STATUS, &optimstatus);
217         check_gurobi_error(grb, error);
218
219         switch (optimstatus) {
220         case GRB_OPTIMAL:           lpp->sol_state = lpp_optimal; break;
221         case GRB_INFEASIBLE:        lpp->sol_state = lpp_infeasible; break;
222         case GRB_INF_OR_UNBD:       lpp->sol_state = lpp_inforunb; break;
223         case GRB_UNBOUNDED:         lpp->sol_state = lpp_unbounded; break;
224         /* TODO: is this correct? */
225         default:                    lpp->sol_state = lpp_feasible; break;
226         }
227
228         if (lpp->sol_state >= lpp_feasible) {
229                 /* get variable solution values */
230                 values = alloca(numcols * sizeof(*values));
231                 error = GRBgetdblattrarray(grb->model, GRB_DBL_ATTR_X, 0, numcols,
232                                            values);
233                 check_gurobi_error(grb, error);
234                 for(i=0; i<numcols; ++i) {
235                         lpp->vars[1+i]->value      = values[i];
236                         lpp->vars[1+i]->value_kind = lpp_value_solution;
237                 }
238
239                 /* Get the value of the objective function. */
240                 error = GRBgetdblattr(grb->model, GRB_DBL_ATTR_OBJVAL, &lpp->objval);
241                 check_gurobi_error(grb, error);
242                 error = GRBgetdblattr(grb->model , GRB_DBL_ATTR_OBJBOUND,
243                                       &lpp->best_bound);
244                 if (error != 0) {
245                         lpp->best_bound = FP_NAN;
246                 }
247         }
248
249         /* get some statistics */
250         error = GRBgetdblattr(grb->model, GRB_DBL_ATTR_ITERCOUNT, &iterations);
251         check_gurobi_error(grb, error);
252         lpp->iterations = (unsigned) iterations;
253
254         error = GRBgetdblattr(grb->model, GRB_DBL_ATTR_RUNTIME, &lpp->sol_time);
255         check_gurobi_error(grb, error);
256 }
257
258 void lpp_solve_gurobi(lpp_t *lpp)
259 {
260         gurobi_t *grb = new_gurobi(lpp);
261         gurobi_construct(grb);
262         gurobi_solve(grb);
263         free_gurobi(grb);
264 }
265
266 #endif