isl_lp.c

   1 /*
   2  * Copyright 2008-2009 Katholieke Universiteit Leuven
   3  *
   4  * Use of this software is governed by the GNU LGPLv2.1 license
   5  *
   6  * Written by Sven Verdoolaege, K.U.Leuven, Departement
   7  * Computerwetenschappen, Celestijnenlaan 200A, B-3001 Leuven, Belgium
   8  */
   9
  10 #include <isl_map_private.h>
  11 #include <isl/ctx.h>
  12 #include <isl/lp.h>
  13 #include "isl_lp_piplib.h"
  14 #include <isl/seq.h>
  15 #include "isl_tab.h"
  16
  17 enum isl_lp_result isl_tab_solve_lp(struct isl_basic_map *bmap, int maximize,
  18                                       isl_int *f, isl_int denom, isl_int *opt,
  19                                       isl_int *opt_denom,
  20                                       struct isl_vec **sol)
  21 {
  22         struct isl_tab *tab;
  23         enum isl_lp_result res;
  24         unsigned dim = isl_basic_map_total_dim(bmap);
  25
  26         if (maximize)
  27                 isl_seq_neg(f, f, 1 + dim);
  28
  29         bmap = isl_basic_map_gauss(bmap, NULL);
  30         tab = isl_tab_from_basic_map(bmap);
  31         res = isl_tab_min(tab, f, denom, opt, opt_denom, 0);
  32         if (res == isl_lp_ok && sol) {
  33                 *sol = isl_tab_get_sample_value(tab);
  34                 if (!*sol)
  35                         res = isl_lp_error;
  36         }
  37         isl_tab_free(tab);
  38
  39         if (maximize)
  40                 isl_seq_neg(f, f, 1 + dim);
  41         if (maximize && opt)
  42                 isl_int_neg(*opt, *opt);
  43
  44         return res;
  45 }
  46
  47 /* Given a basic map "bmap" and an affine combination of the variables "f"
  48  * with denominator "denom", set *opt/*opt_denom to the minimal
  49  * (or maximal if "maximize" is true) value attained by f/d over "bmap",
  50  * assuming the basic map is not empty and the expression cannot attain
  51  * arbitrarily small (or large) values.
  52  * If opt_denom is NULL, then *opt is rounded up (or down)
  53  * to the nearest integer.
  54  * The return value reflects the nature of the result (empty, unbounded,
  55  * minmimal or maximal value returned in *opt).
  56  */
  57 enum isl_lp_result isl_basic_map_solve_lp(struct isl_basic_map *bmap, int max,
  58                                       isl_int *f, isl_int d, isl_int *opt,
  59                                       isl_int *opt_denom,
  60                                       struct isl_vec **sol)
  61 {
  62         if (sol)
  63                 *sol = NULL;
  64
  65         if (!bmap)
  66                 return isl_lp_error;
  67
  68         switch (bmap->ctx->opt->lp_solver) {
  69         case ISL_LP_PIP:
  70                 return isl_pip_solve_lp(bmap, max, f, d, opt, opt_denom, sol);
  71         case ISL_LP_TAB:
  72                 return isl_tab_solve_lp(bmap, max, f, d, opt, opt_denom, sol);
  73         default:
  74                 return isl_lp_error;
  75         }
  76 }
  77
  78 enum isl_lp_result isl_basic_set_solve_lp(struct isl_basic_set *bset, int max,
  79                                       isl_int *f, isl_int d, isl_int *opt,
  80                                       isl_int *opt_denom,
  81                                       struct isl_vec **sol)
  82 {
  83         return isl_basic_map_solve_lp((struct isl_basic_map *)bset, max,
  84                                         f, d, opt, opt_denom, sol);
  85 }
  86
  87 enum isl_lp_result isl_map_solve_lp(__isl_keep isl_map *map, int max,
  88                                       isl_int *f, isl_int d, isl_int *opt,
  89                                       isl_int *opt_denom,
  90                                       struct isl_vec **sol)
  91 {
  92         int i;
  93         isl_int o;
  94         isl_int t;
  95         isl_int opt_i;
  96         isl_int opt_denom_i;
  97         enum isl_lp_result res;
  98         int max_div;
  99         isl_vec *v = NULL;
 100
 101         if (!map)
 102                 return isl_lp_error;
 103         if (map->n == 0)
 104                 return isl_lp_empty;
 105
 106         max_div = 0;
 107         for (i = 0; i < map->n; ++i)
 108                 if (map->p[i]->n_div > max_div)
 109                         max_div = map->p[i]->n_div;
 110         if (max_div > 0) {
 111                 unsigned total = isl_dim_total(map->dim);
 112                 v = isl_vec_alloc(map->ctx, 1 + total + max_div);
 113                 if (!v)
 114                         return isl_lp_error;
 115                 isl_seq_cpy(v->el, f, 1 + total);
 116                 isl_seq_clr(v->el + 1 + total, max_div);
 117                 f = v->el;
 118         }
 119
 120         if (!opt && map->n > 1 && sol) {
 121                 isl_int_init(o);
 122                 opt = &o;
 123         }
 124         if (map->n > 0)
 125                 isl_int_init(opt_i);
 126         if (map->n > 0 && opt_denom) {
 127                 isl_int_init(opt_denom_i);
 128                 isl_int_init(t);
 129         }
 130
 131         res = isl_basic_map_solve_lp(map->p[0], max, f, d,
 132                                         opt, opt_denom, sol);
 133         if (res == isl_lp_error || res == isl_lp_unbounded)
 134                 goto done;
 135
 136         if (sol)
 137                 *sol = NULL;
 138
 139         for (i = 1; i < map->n; ++i) {
 140                 isl_vec *sol_i = NULL;
 141                 enum isl_lp_result res_i;
 142                 int better;
 143
 144                 res_i = isl_basic_map_solve_lp(map->p[i], max, f, d,
 145                                             &opt_i,
 146                                             opt_denom ? &opt_denom_i : NULL,
 147                                             sol ? &sol_i : NULL);
 148                 if (res_i == isl_lp_error || res_i == isl_lp_unbounded) {
 149                         res = res_i;
 150                         goto done;
 151                 }
 152                 if (res_i == isl_lp_empty)
 153                         continue;
 154                 if (res == isl_lp_empty) {
 155                         better = 1;
 156                 } else if (!opt_denom) {
 157                         if (max)
 158                                 better = isl_int_gt(opt_i, *opt);
 159                         else
 160                                 better = isl_int_lt(opt_i, *opt);
 161                 } else {
 162                         isl_int_mul(t, opt_i, *opt_denom);
 163                         isl_int_submul(t, *opt, opt_denom_i);
 164                         if (max)
 165                                 better = isl_int_is_pos(t);
 166                         else
 167                                 better = isl_int_is_neg(t);
 168                 }
 169                 if (better) {
 170                         res = res_i;
 171                         if (opt)
 172                                 isl_int_set(*opt, opt_i);
 173                         if (opt_denom)
 174                                 isl_int_set(*opt_denom, opt_denom_i);
 175                         if (sol) {
 176                                 isl_vec_free(*sol);
 177                                 *sol = sol_i;
 178                         }
 179                 } else
 180                         isl_vec_free(sol_i);
 181         }
 182
 183 done:
 184         isl_vec_free(v);
 185         if (map->n > 0 && opt_denom) {
 186                 isl_int_clear(opt_denom_i);
 187                 isl_int_clear(t);
 188         }
 189         if (map->n > 0)
 190                 isl_int_clear(opt_i);
 191         if (opt == &o)
 192                 isl_int_clear(o);
 193         return res;
 194 }
 195
 196 enum isl_lp_result isl_set_solve_lp(__isl_keep isl_set *set, int max,
 197                                       isl_int *f, isl_int d, isl_int *opt,
 198                                       isl_int *opt_denom,
 199                                       struct isl_vec **sol)
 200 {
 201         return isl_map_solve_lp((struct isl_map *)set, max,
 202                                         f, d, opt, opt_denom, sol);
 203 }