outsourced some functionality
[libfirm] / ir / lpp / lpp.c
1 /**
2  * Author:      Daniel Grund
3  * Date:        Fri 13.05.2005
4  * Copyright:   (c) Universitaet Karlsruhe
5  * Licence:     This file protected by GPL -  GNU GENERAL PUBLIC LICENSE.
6  * CVS-Id:      $Id: lpp.c 27353 2010-04-07 13:33:16Z matze $
7  */
8
9 #include "config.h"
10
11 #ifdef HAVE_IO_H
12 #include <io.h>
13 #endif
14
15 #include <stdlib.h>
16 #include <stdio.h>
17 #include <string.h>
18
19 #include "assert.h"
20 #include "obst.h"
21 #include "hashptr.h"
22 #include "debug.h"
23 #include "set.h"
24
25 #include "sp_matrix.h"
26 #include "mps.h"
27 #include "lpp_t.h"
28 #include "lpp_comm.h"
29 #include "lpp_solvers.h"
30
31 #define HASH_NAME_T(n) HASH_STR((n)->name, strlen((n)->name))
32
33 static firm_dbg_module_t *dbg = NULL;
34
35 static inline char *obst_xstrdup(struct obstack *obst, const char *str) {
36         return obstack_copy0(obst, str, strlen(str));
37 }
38
39 static int cmp_name_t(const void *x, const void *y, size_t size) {
40         const lpp_name_t *n = x;
41         const lpp_name_t *m = y;
42         (void)size; /* stop warnings */
43         return strcmp(n->name, m->name);
44 }
45
46 /**
47  * Update statistic information about matrix usage.
48  */
49 static void update_stats(lpp_t *lpp) {
50         lpp->n_elems    = matrix_get_entries(lpp->m);
51         lpp->matrix_mem = lpp->n_elems * matrix_get_elem_size();
52         lpp->density    = (double)lpp->n_elems / (double)(lpp->cst_next * lpp->var_next) * 100.0;
53 }
54
55 #define INITIAL_SIZE 64
56
57 lpp_t *new_lpp(const char *name, lpp_opt_t opt_type) {
58         return new_lpp_userdef(name, opt_type, INITIAL_SIZE, INITIAL_SIZE, 2.0);
59 }
60
61 lpp_t *new_lpp_userdef(const char *name, lpp_opt_t opt_type,
62                            int estimated_vars, int estimated_csts, double grow_factor)
63 {
64         lpp_t *lpp;
65         int   idx;
66
67         dbg = firm_dbg_register("lpp");
68         lpp = XMALLOCZ(lpp_t);
69         obstack_init(&lpp->obst);
70
71         lpp->name        = obst_xstrdup(&lpp->obst, name);
72         lpp->opt_type    = opt_type;
73         lpp->grow_factor = grow_factor;
74         lpp->cst2nr      = new_set(cmp_name_t, estimated_csts);
75         lpp->var2nr      = new_set(cmp_name_t, estimated_vars);
76         lpp->cst_size    = estimated_csts;
77         lpp->var_size    = estimated_vars;
78         lpp->csts        = XMALLOCNZ(lpp_name_t *, estimated_csts);
79         lpp->vars        = XMALLOCNZ(lpp_name_t *, estimated_vars);
80         lpp->m           = new_matrix(estimated_csts, estimated_vars);
81         lpp->emphasis    = lpp_balanced;
82         idx              = lpp_add_cst(lpp, "obj", lpp_objective, 0);
83         assert(idx == 0);
84         idx              = lpp_add_var(lpp, "rhs", lpp_rhs, 0);
85         assert(idx == 0);
86
87         return lpp;
88 }
89
90 void free_lpp_matrix(lpp_t *lpp) {
91         del_matrix(lpp->m);
92         lpp->m = NULL;
93 }
94
95 void free_lpp(lpp_t *lpp) {
96         obstack_free(&lpp->obst, NULL);
97
98         del_set(lpp->cst2nr);
99         del_set(lpp->var2nr);
100
101         /* matrix might have been already deleted */
102         if (lpp->m)
103                 del_matrix(lpp->m);
104
105         free(lpp->csts);
106         free(lpp->vars);
107
108         free(lpp);
109 }
110
111 double lpp_get_fix_costs(lpp_t *lpp) {
112         return matrix_get(lpp->m, 0, 0);
113 }
114
115 void lpp_set_fix_costs(lpp_t *lpp, double value) {
116         matrix_set(lpp->m, 0, 0, value);
117 }
118
119 static inline int name2nr(set *where, const char *name) {
120         lpp_name_t find, *found;
121         find.name = name;
122         found = set_find(where, &find, sizeof(find), HASH_NAME_T(&find));
123         return (found ? found->nr : -1);
124 }
125
126 #define cst_nr(lpp, name) name2nr(lpp->cst2nr, name)
127 #define var_nr(lpp, name) name2nr(lpp->var2nr, name)
128
129 static inline char *get_next_name(lpp_t *lpp) {
130         char *res = obstack_alloc(&lpp->obst, 12);
131         snprintf(res, 12, "_%d", lpp->next_name_number++);
132         return res;
133 }
134
135 int lpp_add_cst(lpp_t *lpp, const char *cst_name, lpp_cst_t cst_type, double rhs) {
136         lpp_name_t n, *inner;
137
138         DBG((dbg, LEVEL_2, "%s %d %g\n", cst_name, cst_type, rhs));
139
140         if (cst_name && cst_name[0] == '_')
141                 return ERR_NAME_NOT_ALLOWED;
142
143         if (cst_name)
144                 n.name = obst_xstrdup(&lpp->obst, cst_name);
145         else
146                 n.name = get_next_name(lpp);
147
148         n.nr  = -1;
149         inner = set_insert(lpp->cst2nr, &n, sizeof(n), HASH_NAME_T(&n));
150         assert(inner);
151
152         if (inner->nr == -1) {
153                 inner->value_kind    = lpp_none;
154                 inner->value         = 0.0;
155                 inner->nr            = lpp->cst_next;
156                 inner->type.cst_type = cst_type;
157
158                 if (lpp->cst_next == lpp->cst_size) {
159                         lpp->cst_size = (int)((double)lpp->cst_size * lpp->grow_factor) + 1;
160                         lpp->csts     = XREALLOC(lpp->csts, lpp_name_t *, lpp->cst_size);
161                 }
162
163                 lpp->csts[lpp->cst_next] = inner;
164                 lpp->cst_next++;
165                 matrix_set(lpp->m, inner->nr, 0, rhs);
166         }
167
168         update_stats(lpp);
169         return inner->nr;
170 }
171
172 int lpp_add_cst_uniq(lpp_t *lpp, const char *cst_name, lpp_cst_t cst_type, double rhs) {
173         if (cst_name) {
174                 lpp_name_t n;
175
176                 n.name = cst_name;
177                 n.nr   = -1;
178                 assert(!set_find(lpp->cst2nr, &n, sizeof(n), HASH_NAME_T(&n)) &&
179                     "constraint already exists");
180         }
181         return lpp_add_cst(lpp, cst_name, cst_type, rhs);
182 }
183
184 int lpp_get_cst_idx(lpp_t *lpp, const char *cst_name) {
185         DBG((dbg, LEVEL_2, "%s --> %d\n", cst_name, cst_nr(lpp, cst_name)));
186         return cst_nr(lpp, cst_name);
187 }
188
189 void lpp_get_cst_name(lpp_t *lpp, int index, char *buf, size_t buf_size) {
190         DBG((dbg, LEVEL_2, "%d --> %s\n", index, lpp->csts[index]->name));
191         strncpy(buf, lpp->csts[index]->name, buf_size);
192 }
193
194 int lpp_add_var_default(lpp_t *lpp, const char *var_name, lpp_var_t var_type, double obj, double startval) {
195         int val;
196
197         val = lpp_add_var(lpp, var_name, var_type, obj);
198         lpp_set_start_value(lpp, val, startval);
199
200         return val;
201 }
202
203 int lpp_add_var(lpp_t *lpp, const char *var_name, lpp_var_t var_type, double obj) {
204         lpp_name_t n, *inner;
205
206         DBG((dbg, LEVEL_2, "%s %d %g\n", var_name, var_type, obj));
207
208         assert(var_type != lpp_invalid && "invalid is for internal use only");
209
210         if (var_name && var_name[0] == '_')
211                 return ERR_NAME_NOT_ALLOWED;
212
213         if (var_name)
214                 n.name = obst_xstrdup(&lpp->obst, var_name);
215         else
216                 n.name = get_next_name(lpp);
217
218         n.nr  = -1;
219         inner = set_insert(lpp->var2nr, &n, sizeof(n), HASH_NAME_T(&n));
220         assert(inner);
221
222         if (inner->nr == -1) {
223                 inner->nr            = lpp->var_next;
224                 inner->value_kind    = 0;
225                 inner->value         = 0;
226                 inner->type.var_type = var_type;
227
228                 if (lpp->var_next == lpp->var_size) {
229                         lpp->var_size = (int)((double)lpp->var_size * lpp->grow_factor) + 1;
230                         lpp->vars     = XREALLOC(lpp->vars, lpp_name_t *, lpp->var_size);
231                 }
232
233                 lpp->vars[lpp->var_next] = inner;
234                 lpp->var_next++;
235                 matrix_set(lpp->m, 0, inner->nr, obj);
236         }
237
238         update_stats(lpp);
239         return inner->nr;
240 }
241
242 int lpp_get_var_idx(lpp_t *lpp, const char *var_name) {
243         DBG((dbg, LEVEL_2, "%s --> %d\n", var_name, var_nr(lpp, var_name)));
244         return var_nr(lpp, var_name);
245 }
246
247 void lpp_get_var_name(lpp_t *lpp, int index, char *buf, size_t buf_size) {
248         DBG((dbg, LEVEL_2, "%d --> %s\n", index, lpp->vars[index]->name));
249         strncpy(buf, lpp->vars[index]->name, buf_size);
250 }
251
252 int lpp_set_factor(lpp_t *lpp, const char *cst_name, const char *var_name, double value) {
253         int cst, var;
254
255         cst = cst_nr(lpp, cst_name);
256         var = var_nr(lpp, var_name);
257         assert(cst != -1 && var != -1);
258         DBG((dbg, LEVEL_2, "%s[%d] %s[%d] %g\n", cst_name, cst, var_name, var, value));
259         matrix_set(lpp->m, cst, var, value);
260         update_stats(lpp);
261         return 0;
262 }
263
264 int lpp_set_factor_fast(lpp_t *lpp, int cst_idx, int var_idx, double value) {
265         assert(cst_idx >= 0 && var_idx >= 0);
266         assert(cst_idx < lpp->cst_next && var_idx < lpp->var_next);
267         DBG((dbg, LEVEL_2, "%s[%d] %s[%d] %g\n", lpp->csts[cst_idx]->name, cst_idx, lpp->vars[var_idx]->name, var_idx, value));
268         matrix_set(lpp->m, cst_idx, var_idx, value);
269         update_stats(lpp);
270         return 0;
271 }
272
273 int lpp_set_factor_fast_bulk(lpp_t *lpp, int cst_idx, int *var_idx, int num_vars, double value) {
274         assert(cst_idx >= 0 && cst_idx < lpp->cst_next);
275         assert(num_vars < lpp->var_next);
276         DBG((dbg, LEVEL_2, "row %s[%d] %d vars %g\n", lpp->csts[cst_idx]->name, cst_idx, num_vars, value));
277         matrix_set_row_bulk(lpp->m, cst_idx, var_idx, num_vars, value);
278         update_stats(lpp);
279         return 0;
280 }
281
282 void lpp_set_start_value(lpp_t *lpp, int var_idx, double value) {
283         assert(var_idx > 0 && var_idx < lpp->var_next);
284         DBG((dbg, LEVEL_2, "%d %s %g\n", var_idx, lpp->vars[var_idx]->name, value));
285         lpp->vars[var_idx]->value = value;
286         lpp->vars[var_idx]->value_kind = lpp_value_start;
287 }
288
289 lpp_sol_state_t lpp_get_solution(lpp_t *lpp, double *values, int begin, int end) {
290         int i;
291
292         if (lpp->sol_state < lpp_feasible)
293                 return lpp->sol_state;
294
295         /* here we are feasible or optimal */
296         for (i = 0; i < end - begin + 1; ++i)
297                 values[i] = lpp->vars[begin + i]->value;
298
299         return lpp->sol_state;
300 }
301
302 void lpp_check_startvals(lpp_t *lpp) {
303         int cst_idx;
304
305         for (cst_idx = 1; cst_idx < lpp->cst_next; ++cst_idx) {
306                 double     sum     = 0.0;
307                 lpp_name_t *cst    = lpp->csts[cst_idx];
308                 double     cst_val = matrix_get(lpp->m, cst_idx, 0);
309                 int        var_idx;
310
311                 for (var_idx = 1; var_idx < lpp->var_next; ++var_idx) {
312                                 if (lpp->vars[var_idx]->value_kind != lpp_value_start)
313                                         goto next;
314
315                                 sum += lpp->vars[var_idx]->value *
316                                         matrix_get(lpp->m, cst_idx, var_idx);
317                 }
318                 switch(cst->type.cst_type) {
319                         case lpp_equal:
320                                 if(sum != cst_val) {
321                                         fprintf(stderr, "constraint %s unsatisfied: %g != %g\n", cst->name, sum, cst_val);
322                                 }
323                                 break;
324                         case lpp_less:
325                                 if(sum > cst_val) {
326                                         fprintf(stderr, "constraint %s unsatisfied: %g > %g\n", cst->name, sum, cst_val);
327                                 }
328                                 break;
329                         case lpp_greater:
330                                 if(sum < cst_val) {
331                                         fprintf(stderr, "constraint %s unsatisfied: %g < %g\n", cst->name, sum, cst_val);
332                                 }
333                                 break;
334                         default:
335                                 assert(0 && "unknown constraint type");
336                 }
337 next: ;
338         }
339 }
340
341 void lpp_dump(lpp_t *lpp, const char *filename) {
342         FILE *out = fopen(filename, "wt");
343         mps_write_mps(lpp, s_mps_fixed, out);
344         fclose(out);
345 }
346
347 void lpp_set_log(lpp_t *lpp, FILE *log)
348 {
349         lpp->log = log;
350 }
351
352
353 static const char *lpp_cst_op_to_str(lpp_cst_t cst)
354 {
355   switch(cst) {
356     case lpp_equal:
357       return "=";
358     case lpp_less:
359       return "<=";
360     case lpp_greater:
361       return ">=";
362         default:
363           return "";
364   }
365 }
366
367 void lpp_dump_plain(lpp_t *lpp, FILE *f)
368 {
369   int i;
370
371   for(i = 0; i < lpp->cst_next; ++i) {
372     const matrix_elem_t *elm;
373     lpp_name_t *cst = lpp->csts[i];
374
375     fprintf(f, "%16s: ", cst->name);
376     matrix_foreach_in_row(lpp->m, cst->nr, elm) {
377       lpp_name_t *var = lpp->vars[elm->col];
378       /* TODO Perhaps better a define LPP_COL_RHS */
379       if(elm->col > 0)
380         fprintf(f, "%+4.1f*%-16s ", elm->val, var->name);
381     }
382
383     fprintf(f, "%3s %+4.1f\n",
384         lpp_cst_op_to_str(cst->type.cst_type), matrix_get(lpp->m, cst->nr, 0));
385   }
386 }
387
388 /**
389  * Serialize a lpp to a file descriptor.
390  * @param comm The file descriptor.
391  * @param lpp The lpp.
392  */
393 void lpp_serialize(lpp_comm_t *comm, const lpp_t *lpp, int with_names)
394 {
395         int n, i;
396
397         lpp_writel(comm, with_names);
398         lpp_writel(comm, lpp->cst_next);
399         lpp_writel(comm, lpp->var_next);
400         lpp_writel(comm, lpp->opt_type);
401         lpp_writes(comm, lpp->name);
402
403         /* write options */
404         lpp_writel(comm, lpp->set_bound);
405         lpp_writed(comm, lpp->bound);
406         lpp_writed(comm, lpp->time_limit_secs);
407         lpp_writed(comm, lpp->emphasis);
408
409         for(i = 0; i < lpp->cst_next; ++i) {
410                 lpp_name_t *name = lpp->csts[i];
411                 lpp_writel(comm, name->nr);
412                 lpp_writel(comm, name->value_kind);
413                 lpp_writel(comm, name->type.cst_type);
414
415                 if(with_names)
416                         lpp_writes(comm, name->name);
417         }
418
419         for(i = 0; i < lpp->var_next; ++i) {
420                 lpp_name_t *name = lpp->vars[i];
421                 lpp_writel(comm, name->nr);
422                 lpp_writel(comm, name->value_kind);
423                 lpp_writel(comm, name->type.var_type);
424
425                 if(with_names)
426                         lpp_writes(comm, name->name);
427         }
428
429         {
430                 const matrix_elem_t *elm;
431                 n = 0;
432
433                 matrix_foreach(lpp->m, elm)
434                         n++;
435
436                 assert(n == matrix_get_entries(lpp->m));
437                 lpp_writel(comm, n);
438                 matrix_foreach(lpp->m, elm) {
439                         lpp_writel(comm, elm->row);
440                         lpp_writel(comm, elm->col);
441                         lpp_writed(comm, elm->val);
442                 }
443         }
444 }
445
446 #define NAMELEN 64
447
448 /**
449  * Deserialize an lpp from a file descriptor.
450  * @param comm The file descriptor.
451  * @return The Problem.
452  */
453 lpp_t *lpp_deserialize(lpp_comm_t *comm)
454 {
455         int i, n;
456         int with_names;
457
458         lpp_t *lpp = XMALLOCZ(lpp_t);
459
460         /* read general settings */
461         with_names    = lpp_readl(comm);
462         lpp->cst_next = lpp_readl(comm);
463         lpp->var_next = lpp_readl(comm);
464         lpp->opt_type = lpp_readl(comm);
465         lpp->name     = lpp_reads(comm);
466
467         /* read options */
468         lpp->set_bound       = lpp_readl(comm);
469         lpp->bound           = lpp_readd(comm);
470         lpp->time_limit_secs = lpp_readd(comm);
471         lpp->emphasis        = lpp_readd(comm);
472
473         lpp->cst_size = lpp->cst_next;
474         lpp->var_size = lpp->var_next;
475
476         lpp->cst2nr   = new_set(cmp_name_t, lpp->cst_next);
477         lpp->var2nr   = new_set(cmp_name_t, lpp->var_next);
478
479         lpp->csts     = XMALLOCNZ(lpp_name_t*, lpp->cst_next);
480         lpp->vars     = XMALLOCNZ(lpp_name_t*, lpp->var_next);
481         lpp->m        = new_matrix(lpp->cst_next, lpp->var_next);
482
483         for(i = 0; i < lpp->cst_next; ++i) {
484                 lpp_name_t name, *res;
485
486                 name.nr            = lpp_readl(comm);
487                 name.value_kind    = lpp_readl(comm);
488                 name.type.cst_type = lpp_readl(comm);
489
490                 if(with_names)
491                         name.name          = lpp_reads(comm);
492                 else {
493                         char* buf = XMALLOCN(char, NAMELEN);
494                         snprintf(buf, NAMELEN, "c%d\n", name.nr);
495                         name.name = buf;
496                 }
497
498                 res = set_insert(lpp->cst2nr, &name, sizeof(name), HASH_NAME_T(&name));
499                 lpp->csts[name.nr] = res;
500         }
501
502         for(i = 0; i < lpp->var_next; ++i) {
503                 lpp_name_t name, *res;
504
505                 name.nr            = lpp_readl(comm);
506                 name.value_kind    = lpp_readl(comm);
507                 name.type.var_type = lpp_readl(comm);
508
509                 if(with_names)
510                         name.name          = lpp_reads(comm);
511                 else {
512                         char* buf = XMALLOCN(char, NAMELEN);
513                         snprintf(buf, NAMELEN, "v%d\n", name.nr);
514                         name.name = buf;
515                 }
516
517                 res = set_insert(lpp->var2nr, &name, sizeof(name), HASH_NAME_T(&name));
518                 lpp->vars[name.nr] = res;
519         }
520
521         n = lpp_readl(comm);
522         for(i = 0; i < n; ++i) {
523                 matrix_elem_t elm;
524                 elm.row = lpp_readl(comm);
525                 elm.col = lpp_readl(comm);
526                 elm.val = lpp_readd(comm);
527                 matrix_set(lpp->m, elm.row, elm.col, elm.val);
528         }
529
530         return lpp;
531 }
532
533 void lpp_serialize_values(lpp_comm_t *comm, const lpp_t *lpp, lpp_value_kind_t value_kind)
534 {
535   int i, n;
536
537   for(i = 0, n = 0; i < lpp->var_next; ++i)
538     n += lpp->vars[i]->value_kind == value_kind;
539
540   /* Write the number of values to expect */
541   lpp_writel(comm, n);
542
543   /* send the values */
544   for(i = 0, n = lpp->var_next; i < n; ++i) {
545     const lpp_name_t *name = lpp->vars[i];
546     if(name->value_kind == value_kind) {
547       lpp_writel(comm, name->nr);
548       lpp_writed(comm, name->value);
549     }
550   }
551 }
552
553 void lpp_deserialize_values(lpp_comm_t *comm, lpp_t *lpp, lpp_value_kind_t value_kind)
554 {
555   int i, n;
556
557   /* Get the number of values to read */
558   n = lpp_readl(comm);
559
560   for(i = 0; i < n; ++i) {
561     int nr = lpp_readl(comm);
562     lpp_name_t *name = lpp->vars[nr];
563
564     name->value_kind = value_kind;
565     name->value = lpp_readd(comm);
566   }
567 }
568
569 void lpp_serialize_stats(lpp_comm_t *comm, const lpp_t *lpp)
570 {
571         lpp_writel(comm, lpp->sol_state);
572         lpp_writel(comm, lpp->iterations);
573         lpp_writed(comm, lpp->sol_time);
574         lpp_writed(comm, lpp->objval);
575         lpp_writed(comm, lpp->best_bound);
576 }
577
578 void lpp_deserialize_stats(lpp_comm_t *comm, lpp_t *lpp)
579 {
580         lpp->sol_state  = lpp_readl(comm);
581         lpp->iterations = lpp_readl(comm);
582         lpp->sol_time   = lpp_readd(comm);
583         lpp->objval     = lpp_readd(comm);
584         lpp->best_bound = lpp_readd(comm);
585 }
586
587 void lpp_solve(lpp_t *lpp, const char* host, const char* solver)
588 {
589         if (host == NULL || strlen(host) == 0) {
590                 lpp_solver_func_t* f = lpp_find_solver(solver);
591                 if (f != NULL)
592                         f(lpp);
593         }
594
595         else {
596                 lpp_solve_net(lpp, host, solver);
597         }
598 }
599