1 /* 2 * Copyright 2008-2009 Katholieke Universiteit Leuven 3 * Copyright 2010 INRIA Saclay 4 * Copyright 2011 Sven Verdoolaege 5 * 6 * Use of this software is governed by the MIT license 7 * 8 * Written by Sven Verdoolaege, K.U.Leuven, Departement 9 * Computerwetenschappen, Celestijnenlaan 200A, B-3001 Leuven, Belgium 10 * and INRIA Saclay - Ile-de-France, Parc Club Orsay Universite, 11 * ZAC des vignes, 4 rue Jacques Monod, 91893 Orsay, France 12 */ 13 14 #define xSF(TYPE,SUFFIX) TYPE ## SUFFIX 15 #define SF(TYPE,SUFFIX) xSF(TYPE,SUFFIX) 16 17 /* Given a basic map with at least two parallel constraints (as found 18 * by the function parallel_constraints), first look for more constraints 19 * parallel to the two constraint and replace the found list of parallel 20 * constraints by a single constraint with as "input" part the minimum 21 * of the input parts of the list of constraints. Then, recursively call 22 * basic_map_partial_lexopt (possibly finding more parallel constraints) 23 * and plug in the definition of the minimum in the result. 24 * 25 * As in parallel_constraints, only inequality constraints that only 26 * involve input variables that do not occur in any other inequality 27 * constraints are considered. 28 * 29 * More specifically, given a set of constraints 30 * 31 * a x + b_i(p) >= 0 32 * 33 * Replace this set by a single constraint 34 * 35 * a x + u >= 0 36 * 37 * with u a new parameter with constraints 38 * 39 * u <= b_i(p) 40 * 41 * Any solution to the new system is also a solution for the original system 42 * since 43 * 44 * a x >= -u >= -b_i(p) 45 * 46 * Moreover, m = min_i(b_i(p)) satisfies the constraints on u and can 47 * therefore be plugged into the solution. 48 */ 49 static TYPE *SF(basic_map_partial_lexopt_symm,SUFFIX)( 50 __isl_take isl_basic_map *bmap, __isl_take isl_basic_set *dom, 51 __isl_give isl_set **empty, int max, int first, int second) 52 { 53 int i, n, k; 54 int *list = NULL; 55 isl_size bmap_in, bmap_param, bmap_all; 56 unsigned n_in, n_out, n_div; 57 isl_ctx *ctx; 58 isl_vec *var = NULL; 59 isl_mat *cst = NULL; 60 isl_space *map_space, *set_space; 61 62 map_space = isl_basic_map_get_space(bmap); 63 set_space = empty ? isl_basic_set_get_space(dom) : NULL; 64 65 bmap_in = isl_basic_map_dim(bmap, isl_dim_in); 66 bmap_param = isl_basic_map_dim(bmap, isl_dim_param); 67 bmap_all = isl_basic_map_dim(bmap, isl_dim_all); 68 if (bmap_in < 0 || bmap_param < 0 || bmap_all < 0) 69 goto error; 70 n_in = bmap_param + bmap_in; 71 n_out = bmap_all - n_in; 72 73 ctx = isl_basic_map_get_ctx(bmap); 74 list = isl_alloc_array(ctx, int, bmap->n_ineq); 75 var = isl_vec_alloc(ctx, n_out); 76 if ((bmap->n_ineq && !list) || (n_out && !var)) 77 goto error; 78 79 list[0] = first; 80 list[1] = second; 81 isl_seq_cpy(var->el, bmap->ineq[first] + 1 + n_in, n_out); 82 for (i = second + 1, n = 2; i < bmap->n_ineq; ++i) { 83 if (isl_seq_eq(var->el, bmap->ineq[i] + 1 + n_in, n_out) && 84 all_single_occurrence(bmap, i, n_in)) 85 list[n++] = i; 86 } 87 88 cst = isl_mat_alloc(ctx, n, 1 + n_in); 89 if (!cst) 90 goto error; 91 92 for (i = 0; i < n; ++i) 93 isl_seq_cpy(cst->row[i], bmap->ineq[list[i]], 1 + n_in); 94 95 bmap = isl_basic_map_cow(bmap); 96 if (!bmap) 97 goto error; 98 for (i = n - 1; i >= 0; --i) 99 if (isl_basic_map_drop_inequality(bmap, list[i]) < 0) 100 goto error; 101 102 bmap = isl_basic_map_add_dims(bmap, isl_dim_in, 1); 103 bmap = isl_basic_map_extend_constraints(bmap, 0, 1); 104 k = isl_basic_map_alloc_inequality(bmap); 105 if (k < 0) 106 goto error; 107 isl_seq_clr(bmap->ineq[k], 1 + n_in); 108 isl_int_set_si(bmap->ineq[k][1 + n_in], 1); 109 isl_seq_cpy(bmap->ineq[k] + 1 + n_in + 1, var->el, n_out); 110 bmap = isl_basic_map_finalize(bmap); 111 112 n_div = isl_basic_set_dim(dom, isl_dim_div); 113 dom = isl_basic_set_add_dims(dom, isl_dim_set, 1); 114 dom = isl_basic_set_extend_constraints(dom, 0, n); 115 for (i = 0; i < n; ++i) { 116 k = isl_basic_set_alloc_inequality(dom); 117 if (k < 0) 118 goto error; 119 isl_seq_cpy(dom->ineq[k], cst->row[i], 1 + n_in); 120 isl_int_set_si(dom->ineq[k][1 + n_in], -1); 121 isl_seq_clr(dom->ineq[k] + 1 + n_in + 1, n_div); 122 } 123 124 isl_vec_free(var); 125 free(list); 126 127 return SF(basic_map_partial_lexopt_symm_core,SUFFIX)(bmap, dom, empty, 128 max, cst, map_space, set_space); 129 error: 130 isl_space_free(map_space); 131 isl_space_free(set_space); 132 isl_mat_free(cst); 133 isl_vec_free(var); 134 free(list); 135 isl_basic_set_free(dom); 136 isl_basic_map_free(bmap); 137 return NULL; 138 } 139 140 /* Recursive part of isl_tab_basic_map_partial_lexopt*, after detecting 141 * equalities and removing redundant constraints. 142 * 143 * We first check if there are any parallel constraints (left). 144 * If not, we are in the base case. 145 * If there are parallel constraints, we replace them by a single 146 * constraint in basic_map_partial_lexopt_symm_pma and then call 147 * this function recursively to look for more parallel constraints. 148 */ 149 static __isl_give TYPE *SF(basic_map_partial_lexopt,SUFFIX)( 150 __isl_take isl_basic_map *bmap, __isl_take isl_basic_set *dom, 151 __isl_give isl_set **empty, int max) 152 { 153 isl_bool par = isl_bool_false; 154 int first, second; 155 156 if (!bmap) 157 goto error; 158 159 if (bmap->ctx->opt->pip_symmetry) 160 par = parallel_constraints(bmap, &first, &second); 161 if (par < 0) 162 goto error; 163 if (!par) 164 return SF(basic_map_partial_lexopt_base,SUFFIX)(bmap, dom, 165 empty, max); 166 167 return SF(basic_map_partial_lexopt_symm,SUFFIX)(bmap, dom, empty, max, 168 first, second); 169 error: 170 isl_basic_set_free(dom); 171 isl_basic_map_free(bmap); 172 return NULL; 173 } 174 175 /* Compute the lexicographic minimum (or maximum if "flags" includes 176 * ISL_OPT_MAX) of "bmap" over the domain "dom" and return the result as 177 * either a map or a piecewise multi-affine expression depending on TYPE. 178 * If "empty" is not NULL, then *empty is assigned a set that 179 * contains those parts of the domain where there is no solution. 180 * If "flags" includes ISL_OPT_FULL, then "dom" is NULL and the optimum 181 * should be computed over the domain of "bmap". "empty" is also NULL 182 * in this case. 183 * If "bmap" is marked as rational (ISL_BASIC_MAP_RATIONAL), 184 * then we compute the rational optimum. Otherwise, we compute 185 * the integral optimum. 186 * 187 * We perform some preprocessing. As the PILP solver does not 188 * handle implicit equalities very well, we first make sure all 189 * the equalities are explicitly available. 190 * 191 * We also add context constraints to the basic map and remove 192 * redundant constraints. This is only needed because of the 193 * way we handle simple symmetries. In particular, we currently look 194 * for symmetries on the constraints, before we set up the main tableau. 195 * It is then no good to look for symmetries on possibly redundant constraints. 196 * If the domain was extracted from the basic map, then there is 197 * no need to add back those constraints again. 198 */ 199 __isl_give TYPE *SF(isl_tab_basic_map_partial_lexopt,SUFFIX)( 200 __isl_take isl_basic_map *bmap, __isl_take isl_basic_set *dom, 201 __isl_give isl_set **empty, unsigned flags) 202 { 203 int max, full; 204 isl_bool compatible; 205 206 if (empty) 207 *empty = NULL; 208 209 full = ISL_FL_ISSET(flags, ISL_OPT_FULL); 210 if (full) 211 dom = extract_domain(bmap, flags); 212 compatible = isl_basic_map_compatible_domain(bmap, dom); 213 if (compatible < 0) 214 goto error; 215 if (!compatible) 216 isl_die(isl_basic_map_get_ctx(bmap), isl_error_invalid, 217 "domain does not match input", goto error); 218 219 max = ISL_FL_ISSET(flags, ISL_OPT_MAX); 220 if (isl_basic_set_dim(dom, isl_dim_all) == 0) 221 return SF(basic_map_partial_lexopt,SUFFIX)(bmap, dom, empty, 222 max); 223 224 if (!full) 225 bmap = isl_basic_map_intersect_domain(bmap, 226 isl_basic_set_copy(dom)); 227 bmap = isl_basic_map_detect_equalities(bmap); 228 bmap = isl_basic_map_remove_redundancies(bmap); 229 230 return SF(basic_map_partial_lexopt,SUFFIX)(bmap, dom, empty, max); 231 error: 232 isl_basic_set_free(dom); 233 isl_basic_map_free(bmap); 234 return NULL; 235 } 236