9 struct isl_tab_undo **stack;
17 /* number of constraints in initial product tableau */
19 /* objective function has fixed or no integer value */
23 static struct tab_lp *init_lp(struct isl_tab *tab);
24 static void set_lp_obj(struct tab_lp *lp, isl_int *row, int dim);
25 static int solve_lp(struct tab_lp *lp);
26 static void get_obj_val(struct tab_lp* lp, mpq_t *F);
27 static void delete_lp(struct tab_lp *lp);
28 static int add_lp_row(struct tab_lp *lp, isl_int *row, int dim);
29 static void get_alpha(struct tab_lp* lp, int row, mpq_t *alpha);
30 static void del_lp_row(struct tab_lp *lp);
31 static int cut_lp_to_hyperplane(struct tab_lp *lp, isl_int *row);
33 #define GBR_LP struct tab_lp
34 #define GBR_type mpq_t
35 #define GBR_init(v) mpq_init(v)
36 #define GBR_clear(v) mpq_clear(v)
37 #define GBR_set(a,b) mpq_set(a,b)
38 #define GBR_set_ui(a,b) mpq_set_ui(a,b,1)
39 #define GBR_mul(a,b,c) mpq_mul(a,b,c)
40 #define GBR_lt(a,b) (mpq_cmp(a,b) < 0)
41 #define GBR_is_zero(a) (mpq_sgn(a) == 0)
42 #define GBR_floor(a,b) mpz_fdiv_q(a,mpq_numref(b),mpq_denref(b))
43 #define GBR_ceil(a,b) mpz_cdiv_q(a,mpq_numref(b),mpq_denref(b))
44 #define GBR_lp_init(P) init_lp(P)
45 #define GBR_lp_set_obj(lp, obj, dim) set_lp_obj(lp, obj, dim)
46 #define GBR_lp_solve(lp) solve_lp(lp)
47 #define GBR_lp_get_obj_val(lp, F) get_obj_val(lp, F)
48 #define GBR_lp_delete(lp) delete_lp(lp)
49 #define GBR_lp_next_row(lp) lp->neq
50 #define GBR_lp_add_row(lp, row, dim) add_lp_row(lp, row, dim)
51 #define GBR_lp_get_alpha(lp, row, alpha) get_alpha(lp, row, alpha)
52 #define GBR_lp_del_row(lp) del_lp_row(lp)
53 #define GBR_lp_is_fixed(lp) (lp)->is_fixed
54 #define GBR_lp_cut(lp, obj) cut_lp_to_hyperplane(lp, obj)
55 #include "basis_reduction_templ.c"
57 /* Set up a tableau for the Cartesian product of bset with itself.
58 * This could be optimized by first setting up a tableau for bset
59 * and then performing the Cartesian product on the tableau.
61 static struct isl_tab *gbr_tab(struct isl_tab *tab, struct isl_vec *row)
71 prod = isl_tab_product(tab, tab);
72 if (isl_tab_extend_cons(prod, 3 * dim + 1) < 0) {
79 static struct tab_lp *init_lp(struct isl_tab *tab)
81 struct tab_lp *lp = NULL;
86 lp = isl_calloc_type(tab->mat->ctx, struct tab_lp);
90 isl_int_init(lp->opt);
91 isl_int_init(lp->opt_denom);
92 isl_int_init(lp->tmp);
93 isl_int_init(lp->tmp2);
97 lp->ctx = tab->mat->ctx;
100 lp->stack = isl_alloc_array(lp->ctx, struct isl_tab_undo *, lp->dim);
102 lp->row = isl_vec_alloc(lp->ctx, 1 + 2 * lp->dim);
105 lp->tab = gbr_tab(tab, lp->row);
108 lp->con_offset = lp->tab->n_con;
118 static void set_lp_obj(struct tab_lp *lp, isl_int *row, int dim)
123 static int solve_lp(struct tab_lp *lp)
125 enum isl_lp_result res;
130 isl_int_set_si(lp->row->el[0], 0);
131 isl_seq_cpy(lp->row->el + 1, lp->obj, lp->dim);
132 isl_seq_neg(lp->row->el + 1 + lp->dim, lp->obj, lp->dim);
134 flags = ISL_TAB_SAVE_DUAL;
135 res = isl_tab_min(lp->tab, lp->row->el, lp->ctx->one,
136 &lp->opt, &lp->opt_denom, flags);
137 isl_int_mul_ui(lp->opt_denom, lp->opt_denom, 2);
138 if (isl_int_abs_lt(lp->opt, lp->opt_denom)) {
139 struct isl_vec *sample = isl_tab_get_sample_value(lp->tab);
142 isl_seq_inner_product(lp->obj, sample->el + 1, lp->dim, &lp->tmp);
143 isl_seq_inner_product(lp->obj, sample->el + 1 + lp->dim, lp->dim, &lp->tmp2);
144 isl_int_cdiv_q(lp->tmp, lp->tmp, sample->el[0]);
145 isl_int_fdiv_q(lp->tmp2, lp->tmp2, sample->el[0]);
146 if (isl_int_ge(lp->tmp, lp->tmp2))
148 isl_vec_free(sample);
150 isl_int_divexact_ui(lp->opt_denom, lp->opt_denom, 2);
151 if (res != isl_lp_ok)
156 /* The current objective function has a fixed (or no) integer value.
157 * Cut the tableau to the hyperplane that fixes this value in
158 * both halves of the tableau.
159 * Return 1 if the resulting tableau is empty.
161 static int cut_lp_to_hyperplane(struct tab_lp *lp, isl_int *row)
163 enum isl_lp_result res;
165 isl_int_set_si(lp->row->el[0], 0);
166 isl_seq_cpy(lp->row->el + 1, row, lp->dim);
167 isl_seq_clr(lp->row->el + 1 + lp->dim, lp->dim);
168 res = isl_tab_min(lp->tab, lp->row->el, lp->ctx->one,
170 if (res != isl_lp_ok)
173 isl_int_neg(lp->row->el[0], lp->tmp);
174 lp->tab = isl_tab_add_eq(lp->tab, lp->row->el);
176 isl_seq_cpy(lp->row->el + 1 + lp->dim, row, lp->dim);
177 isl_seq_clr(lp->row->el + 1, lp->dim);
178 lp->tab = isl_tab_add_eq(lp->tab, lp->row->el);
182 return lp->tab->empty;
185 static void get_obj_val(struct tab_lp* lp, mpq_t *F)
187 isl_int_neg(mpq_numref(*F), lp->opt);
188 isl_int_set(mpq_denref(*F), lp->opt_denom);
191 static void delete_lp(struct tab_lp *lp)
196 isl_int_clear(lp->opt);
197 isl_int_clear(lp->opt_denom);
198 isl_int_clear(lp->tmp);
199 isl_int_clear(lp->tmp2);
200 isl_vec_free(lp->row);
202 isl_tab_free(lp->tab);
203 isl_ctx_deref(lp->ctx);
207 static int add_lp_row(struct tab_lp *lp, isl_int *row, int dim)
209 lp->stack[lp->neq] = isl_tab_snap(lp->tab);
211 isl_int_set_si(lp->row->el[0], 0);
212 isl_seq_cpy(lp->row->el + 1, row, lp->dim);
213 isl_seq_neg(lp->row->el + 1 + lp->dim, row, lp->dim);
215 lp->tab = isl_tab_add_valid_eq(lp->tab, lp->row->el);
220 static void get_alpha(struct tab_lp* lp, int row, mpq_t *alpha)
222 row += lp->con_offset;
223 isl_int_neg(mpq_numref(*alpha), lp->tab->dual->el[1 + row]);
224 isl_int_set(mpq_denref(*alpha), lp->tab->dual->el[0]);
227 static void del_lp_row(struct tab_lp *lp)
230 isl_tab_rollback(lp->tab, lp->stack[lp->neq]);