consistently give all lpp functions an lpp_ prefix
[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
32 #include "obst.h"
33
34 #include <gurobi_c.h>
35
36 #include "error.h"
37 #include "sp_matrix.h"
38
39 static char gurobi_cst_encoding[4] = { 0, GRB_EQUAL, GRB_LESS_EQUAL, GRB_GREATER_EQUAL };
40 static char gurobi_var_encoding[4] = { 0, 0, GRB_CONTINUOUS, GRB_BINARY };
41
42 typedef struct _gurobi_t {
43         lpp_t *lpp;
44         GRBenv *env;
45         GRBmodel *model;
46 } gurobi_t;
47
48 static void check_gurobi_error(gurobi_t *grb, int error)
49 {
50         if (error != 0) {
51                 panic("gurobi error: %s", GRBgeterrormsg(grb->env));
52         }
53 }
54
55 static gurobi_t *new_gurobi(lpp_t *lpp)
56 {
57         int error;
58
59         gurobi_t *grb = XMALLOCZ(gurobi_t);
60         grb->lpp = lpp;
61         error = GRBloadenv(&grb->env, NULL);
62         check_gurobi_error(grb, error);
63         error = GRBsetlogfile(grb->env, lpp->log);
64         check_gurobi_error(grb, error);
65
66         return grb;
67 }
68
69 static void free_gurobi(gurobi_t *grb)
70 {
71         GRBfreeenv(grb->env);
72         free(grb);
73 }
74
75 /**
76  * Build CPLEX data structure from LPP matrix.
77  * @note: The LPP matrix is freed after this step, to save memory.
78  */
79 static void gurobi_construct(gurobi_t *grb)
80 {
81         const matrix_elem_t *elem;
82         int                  i, o;
83         //int                  sv_cnt;
84         //int                 *indices;
85         //double              *startv;
86         int                  numcols, numrows, numentries;
87         int                  objsen, *matbeg, *matcnt, *matind;
88         double               *obj, *rhs, *matval, *lb;
89         char                 *sense, *vartype;
90         char                 **colname, **rowname;
91         struct obstack       obst;
92         lpp_t                *lpp = grb->lpp;
93         int                  error;
94
95         numcols    = lpp->var_next-1;
96         numrows    = lpp->cst_next-1;
97         numentries = matrix_get_entries(lpp->m);
98         objsen     = lpp->opt_type == lpp_minimize ? 1 : -1;
99         obstack_init(&obst);
100
101         obj     = obstack_alloc(&obst, numcols * sizeof(*obj));
102         lb      = obstack_alloc(&obst, numcols * sizeof(*lb));
103         colname = obstack_alloc(&obst, numcols * sizeof(*colname));
104         rowname = obstack_alloc(&obst, numrows * sizeof(*rowname));
105         vartype = obstack_alloc(&obst, numcols * sizeof(*vartype));
106         //indices = obstack_alloc(&obst, numcols * sizeof(*indices));
107         //startv  = obstack_alloc(&obst, numcols * sizeof(*startv));
108         matbeg  = obstack_alloc(&obst, numcols * sizeof(*matbeg));
109         matcnt  = obstack_alloc(&obst, numcols * sizeof(*matcnt));
110         matind  = obstack_alloc(&obst, numentries * sizeof(*matind));
111         matval  = obstack_alloc(&obst, numentries * sizeof(*matval));
112         rhs     = obstack_alloc(&obst, numrows * sizeof(*rhs));
113         sense   = obstack_alloc(&obst, numrows * sizeof(*sense));
114
115         o      = 0;
116         //sv_cnt = 0;
117         /* fill the CPLEX matrix*/
118         for (i = 0; i < numcols; ++i) {
119                 lpp_name_t *curr_var = lpp->vars[1+i];
120
121                 obj[i] = matrix_get(lpp->m, 0, 1+i);
122                 lb[i]  = 0.0;
123
124                 colname[i] = (char*) curr_var->name;
125                 vartype[i] = gurobi_var_encoding[curr_var->type.var_type];
126
127 #if 0
128                 if (curr_var->value_kind == lpp_value_start) {
129                         panic("start values not supported in gurobi yet");
130                         indices[sv_cnt]  = i;
131                         startv[sv_cnt++] = curr_var->value;
132                 }
133 #endif
134
135                 matbeg[i] = o;
136                 matcnt[i] = 0;
137                 matrix_foreach_in_col(lpp->m, 1 + i, elem) {
138                         if (elem->row == 0)
139                                 continue;
140                         matind[o] = elem->row-1;
141                         matval[o] = elem->val;
142                         matcnt[i]++;
143                         o++;
144                 }
145         }
146
147         /* get constraint stuff (right hand side, type, name) */
148         for (i = 0; i < numrows; ++i) {
149                 lpp_name_t *curr_cst = lpp->csts[1 + i];
150
151                 rhs[i]     = matrix_get(lpp->m, 1 + i, 0);
152                 sense[i]   = gurobi_cst_encoding[curr_cst->type.cst_type];
153                 rowname[i] = (char*) curr_cst->name;
154         }
155
156         error = GRBloadmodel(grb->env, &grb->model, lpp->name, numcols, numrows,
157                              objsen, 0, obj, sense, rhs, matbeg, matcnt, matind,
158                              matval, lb, NULL, vartype, colname, rowname);
159         check_gurobi_error(grb, error);
160
161         obstack_free(&obst, NULL);
162         lpp_free_matrix(lpp);
163 }
164
165 static void gurobi_solve(gurobi_t *grb)
166 {
167         lpp_t *lpp = grb->lpp;
168         int i;
169         int optimstatus;
170         int error;
171         int numcols = lpp->var_next-1;
172         double *values;
173         double  iterations;
174
175         /* set performance parameters */
176         // CPXsetintparam(grb->env, CPX_PARAM_MIPSTART, CPX_ON);
177         //CPXsetintparam(grb->env, CPX_PARAM_MIPORDTYPE, CPX_MIPORDER_COST);
178         /* output every search tree node */
179         // CPXsetintparam(grb->env, CPX_PARAM_MIPINTERVAL, 1);
180
181         /* experimental switches */
182         // CPXsetintparam(grb->env, CPX_PARAM_VARSEL, CPX_VARSEL_STRONG);
183         // CPXsetdblparam(grb->env, CPX_PARAM_BTTOL, 1.0);
184         // CPXsetintparam(grb->env, CPX_PARAM_BRDIR, CPX_BRDIR_UP);
185
186         /* Set the time limit appropriately */
187         if(lpp->time_limit_secs > 0.0) {
188                 error = GRBsetdblparam(grb->env, GRB_DBL_PAR_TIMELIMIT, lpp->time_limit_secs);
189                 check_gurobi_error(grb, error);
190         }
191
192         /*
193          * If we have enough time, we instruct cplex to imply some
194          * of its higher order magic to pursue the best solution
195          */
196         if(lpp->emphasis) {
197                 /* not implemented */
198         }
199
200         /*
201          * If a bound of the objective function is supplied,
202          * set it accordingly, dependign on minimization or maximization.
203          */
204         if(lpp->set_bound) {
205                 //panic("bound not implemented yet");
206                 fprintf(stderr, "Warning: gurobi bound not implemented yet\n");
207         }
208
209         /* solve */
210         error = GRBoptimize(grb->model);
211         check_gurobi_error(grb, error);
212
213         /* get solution status */
214         error = GRBgetintattr(grb->model, GRB_INT_ATTR_STATUS, &optimstatus);
215         check_gurobi_error(grb, error);
216
217         switch (optimstatus) {
218         case GRB_OPTIMAL:           lpp->sol_state = lpp_optimal; break;
219         case GRB_INFEASIBLE:        lpp->sol_state = lpp_infeasible; break;
220         case GRB_INF_OR_UNBD:       lpp->sol_state = lpp_inforunb; break;
221         case GRB_UNBOUNDED:         lpp->sol_state = lpp_unbounded; break;
222         /* TODO: is this correct? */
223         default:                    lpp->sol_state = lpp_feasible; break;
224         }
225
226         /* get variable solution values */
227         values = alloca(numcols * sizeof(*values));
228         error = GRBgetdblattrarray(grb->model, GRB_DBL_ATTR_X, 0, numcols, values);
229         check_gurobi_error(grb, error);
230         for(i=0; i<numcols; ++i) {
231                 lpp->vars[1+i]->value      = values[i];
232                 lpp->vars[1+i]->value_kind = lpp_value_solution;
233         }
234
235         /* Get the value of the objective function. */
236         error = GRBgetdblattr(grb->model, GRB_DBL_ATTR_OBJVAL, &lpp->objval);
237         check_gurobi_error(grb, error);
238         error = GRBgetdblattr(grb->model , GRB_DBL_ATTR_OBJBOUND, &lpp->best_bound);
239         check_gurobi_error(grb, error);
240
241         /* get some statistics */
242         error = GRBgetdblattr(grb->model, GRB_DBL_ATTR_ITERCOUNT, &iterations);
243         check_gurobi_error(grb, error);
244         lpp->iterations = (unsigned) iterations;
245
246         error = GRBgetdblattr(grb->model, GRB_DBL_ATTR_RUNTIME, &lpp->sol_time);
247         check_gurobi_error(grb, error);
248 }
249
250 void lpp_solve_gurobi(lpp_t *lpp)
251 {
252         gurobi_t *grb = new_gurobi(lpp);
253         gurobi_construct(grb);
254         gurobi_solve(grb);
255         free_gurobi(grb);
256 }
257
258 #endif