isl_pw_qpolynomial_fold_bound: simplify construction of target space
[platform/upstream/isl.git] / pip.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 <assert.h>
11 #include <string.h>
12 #include <strings.h>
13 #include <isl_map_private.h>
14 #include <isl/set.h>
15 #include "isl_tab.h"
16 #include "isl_sample.h"
17 #include "isl_scan.h"
18 #include <isl/seq.h>
19 #include <isl/ilp.h>
20 #include <isl_point_private.h>
21
22 /* The input of this program is the same as that of the "example" program
23  * from the PipLib distribution, except that the "big parameter column"
24  * should always be -1.
25  *
26  * Context constraints in PolyLib format
27  * -1
28  * Problem constraints in PolyLib format
29  * Optional list of options
30  *
31  * The options are
32  *      Maximize        compute maximum instead of minimum
33  *      Rational        compute rational optimum instead of integer optimum
34  *      Urs_parms       don't assume parameters are non-negative
35  *      Urs_unknowns    don't assume unknowns are non-negative
36  */
37
38 struct options {
39         struct isl_options      *isl;
40         unsigned                 verify;
41 };
42
43 struct isl_arg options_arg[] = {
44 ISL_ARG_CHILD(struct options, isl, "isl", isl_options_arg, "isl options")
45 ISL_ARG_BOOL(struct options, verify, 'T', "verify", 0, NULL)
46 ISL_ARG_END
47 };
48
49 ISL_ARG_DEF(options, struct options, options_arg)
50
51 static __isl_give isl_basic_set *set_bounds(__isl_take isl_basic_set *bset)
52 {
53         unsigned nparam;
54         int i, r;
55         isl_point *pt, *pt2;
56         isl_basic_set *box;
57
58         nparam = isl_basic_set_dim(bset, isl_dim_param);
59         r = nparam >= 8 ? 4 : nparam >= 5 ? 6 : 30;
60
61         pt = isl_basic_set_sample_point(isl_basic_set_copy(bset));
62         pt2 = isl_point_copy(pt);
63
64         for (i = 0; i < nparam; ++i) {
65                 pt = isl_point_add_ui(pt, isl_dim_param, i, r);
66                 pt2 = isl_point_sub_ui(pt2, isl_dim_param, i, r);
67         }
68
69         box = isl_basic_set_box_from_points(pt, pt2);
70
71         return isl_basic_set_intersect(bset, box);
72 }
73
74 static struct isl_basic_set *to_parameter_domain(struct isl_basic_set *context)
75 {
76         return isl_basic_set_move_dims(context, isl_dim_param, 0, isl_dim_set, 0,
77                                        isl_basic_set_dim(context, isl_dim_set));
78 }
79
80 isl_basic_set *plug_in_parameters(isl_basic_set *bset, struct isl_vec *params)
81 {
82         int i;
83
84         for (i = 0; i < params->size - 1; ++i)
85                 bset = isl_basic_set_fix(bset,
86                                          isl_dim_param, i, params->el[1 + i]);
87
88         bset = isl_basic_set_remove_dims(bset,
89                                          isl_dim_param, 0, params->size - 1);
90
91         isl_vec_free(params);
92
93         return bset;
94 }
95
96 isl_set *set_plug_in_parameters(isl_set *set, struct isl_vec *params)
97 {
98         int i;
99
100         for (i = 0; i < params->size - 1; ++i)
101                 set = isl_set_fix(set, isl_dim_param, i, params->el[1 + i]);
102
103         set = isl_set_remove_dims(set, isl_dim_param, 0, params->size - 1);
104
105         isl_vec_free(params);
106
107         return set;
108 }
109
110 /* Compute the lexicographically minimal (or maximal if max is set)
111  * element of bset for the given values of the parameters, by
112  * successively solving an ilp problem in each direction.
113  */
114 struct isl_vec *opt_at(struct isl_basic_set *bset,
115         struct isl_vec *params, int max)
116 {
117         unsigned dim;
118         struct isl_vec *opt;
119         struct isl_vec *obj;
120         int i;
121
122         dim = isl_basic_set_dim(bset, isl_dim_set);
123
124         bset = plug_in_parameters(bset, params);
125
126         if (isl_basic_set_plain_is_empty(bset)) {
127                 opt = isl_vec_alloc(bset->ctx, 0);
128                 isl_basic_set_free(bset);
129                 return opt;
130         }
131
132         opt = isl_vec_alloc(bset->ctx, 1 + dim);
133         assert(opt);
134
135         obj = isl_vec_alloc(bset->ctx, 1 + dim);
136         assert(obj);
137
138         isl_int_set_si(opt->el[0], 1);
139         isl_int_set_si(obj->el[0], 0);
140
141         for (i = 0; i < dim; ++i) {
142                 enum isl_lp_result res;
143
144                 isl_seq_clr(obj->el + 1, dim);
145                 isl_int_set_si(obj->el[1 + i], 1);
146                 res = isl_basic_set_solve_ilp(bset, max, obj->el,
147                                                 &opt->el[1 + i], NULL);
148                 if (res == isl_lp_empty)
149                         goto empty;
150                 assert(res == isl_lp_ok);
151                 bset = isl_basic_set_fix(bset, isl_dim_set, i, opt->el[1 + i]);
152         }
153
154         isl_basic_set_free(bset);
155         isl_vec_free(obj);
156
157         return opt;
158 empty:
159         isl_vec_free(opt);
160         opt = isl_vec_alloc(bset->ctx, 0);
161         isl_basic_set_free(bset);
162         isl_vec_free(obj);
163
164         return opt;
165 }
166
167 struct isl_scan_pip {
168         struct isl_scan_callback callback;
169         isl_basic_set *bset;
170         isl_set *sol;
171         isl_set *empty;
172         int stride;
173         int n;
174         int max;
175 };
176
177 /* Check if the "manually" computed optimum of bset at the "sample"
178  * values of the parameters agrees with the solution of pilp problem
179  * represented by the pair (sol, empty).
180  * In particular, if there is no solution for this value of the parameters,
181  * then it should be an element of the parameter domain "empty".
182  * Otherwise, the optimal solution, should be equal to the result of
183  * plugging in the value of the parameters in "sol".
184  */
185 static int scan_one(struct isl_scan_callback *callback,
186         __isl_take isl_vec *sample)
187 {
188         struct isl_scan_pip *sp = (struct isl_scan_pip *)callback;
189         struct isl_vec *opt;
190
191         sp->n--;
192
193         opt = opt_at(isl_basic_set_copy(sp->bset), isl_vec_copy(sample), sp->max);
194         assert(opt);
195
196         if (opt->size == 0) {
197                 isl_point *sample_pnt;
198                 sample_pnt = isl_point_alloc(isl_set_get_space(sp->empty), sample);
199                 assert(isl_set_contains_point(sp->empty, sample_pnt));
200                 isl_point_free(sample_pnt);
201                 isl_vec_free(opt);
202         } else {
203                 isl_set *sol;
204                 isl_set *opt_set;
205                 opt_set = isl_set_from_basic_set(isl_basic_set_from_vec(opt));
206                 sol = set_plug_in_parameters(isl_set_copy(sp->sol), sample);
207                 assert(isl_set_is_equal(opt_set, sol));
208                 isl_set_free(sol);
209                 isl_set_free(opt_set);
210         }
211
212         if (!(sp->n % sp->stride)) {
213                 printf("o");
214                 fflush(stdout);
215         }
216
217         return sp->n >= 1 ? 0 : -1;
218 }
219
220 static void check_solution(isl_basic_set *bset, isl_basic_set *context,
221         isl_set *sol, isl_set *empty, int max)
222 {
223         struct isl_scan_pip sp;
224         isl_int count, count_max;
225         int i, n;
226         int r;
227
228         context = set_bounds(context);
229         context = isl_basic_set_underlying_set(context);
230
231         isl_int_init(count);
232         isl_int_init(count_max);
233
234         isl_int_set_si(count_max, 2000);
235         r = isl_basic_set_count_upto(context, count_max, &count);
236         assert(r >= 0);
237         n = isl_int_get_si(count);
238
239         isl_int_clear(count_max);
240         isl_int_clear(count);
241
242         sp.callback.add = scan_one;
243         sp.bset = bset;
244         sp.sol = sol;
245         sp.empty = empty;
246         sp.n = n;
247         sp.stride = n > 70 ? 1 + (n + 1)/70 : 1;
248         sp.max = max;
249
250         for (i = 0; i < n; i += sp.stride)
251                 printf(".");
252         printf("\r");
253         fflush(stdout);
254
255         isl_basic_set_scan(context, &sp.callback);
256
257         printf("\n");
258
259         isl_basic_set_free(bset);
260 }
261
262 int main(int argc, char **argv)
263 {
264         struct isl_ctx *ctx;
265         struct isl_basic_set *context, *bset, *copy, *context_copy;
266         struct isl_set *set;
267         struct isl_set *empty;
268         int neg_one;
269         char s[1024];
270         int urs_parms = 0;
271         int urs_unknowns = 0;
272         int max = 0;
273         int rational = 0;
274         int n;
275         struct options *options;
276
277         options = options_new_with_defaults();
278         assert(options);
279         argc = options_parse(options, argc, argv, ISL_ARG_ALL);
280
281         ctx = isl_ctx_alloc_with_options(options_arg, options);
282
283         context = isl_basic_set_read_from_file(ctx, stdin, 0);
284         assert(context);
285         n = fscanf(stdin, "%d", &neg_one);
286         assert(n == 1);
287         assert(neg_one == -1);
288         bset = isl_basic_set_read_from_file(ctx, stdin,
289                 isl_basic_set_dim(context, isl_dim_set));
290
291         while (fgets(s, sizeof(s), stdin)) {
292                 if (strncasecmp(s, "Maximize", 8) == 0)
293                         max = 1;
294                 if (strncasecmp(s, "Rational", 8) == 0) {
295                         rational = 1;
296                         bset = isl_basic_set_set_rational(bset);
297                 }
298                 if (strncasecmp(s, "Urs_parms", 9) == 0)
299                         urs_parms = 1;
300                 if (strncasecmp(s, "Urs_unknowns", 12) == 0)
301                         urs_unknowns = 1;
302         }
303         if (!urs_parms)
304                 context = isl_basic_set_intersect(context,
305                 isl_basic_set_positive_orthant(isl_basic_set_get_space(context)));
306         context = to_parameter_domain(context);
307         if (!urs_unknowns)
308                 bset = isl_basic_set_intersect(bset,
309                 isl_basic_set_positive_orthant(isl_basic_set_get_space(bset)));
310
311         if (options->verify) {
312                 copy = isl_basic_set_copy(bset);
313                 context_copy = isl_basic_set_copy(context);
314         }
315
316         if (max)
317                 set = isl_basic_set_partial_lexmax(bset, context, &empty);
318         else
319                 set = isl_basic_set_partial_lexmin(bset, context, &empty);
320
321         if (options->verify) {
322                 assert(!rational);
323                 check_solution(copy, context_copy, set, empty, max);
324         } else {
325                 isl_set_print(set, stdout, 0, ISL_FORMAT_ISL);
326                 fprintf(stdout, "\n");
327                 fprintf(stdout, "no solution: ");
328                 isl_set_print(empty, stdout, 0, ISL_FORMAT_ISL);
329                 fprintf(stdout, "\n");
330         }
331
332         isl_set_free(set);
333         isl_set_free(empty);
334         isl_ctx_free(ctx);
335
336         return 0;
337 }