isl_map.c: parameter_compute_divs: check result of isl_basic_set_gauss call
[isl.git] / pip.c
blobda08d2e19ac97a1d69d1e331de680610a425cc2a
1 /*
2 * Copyright 2008-2009 Katholieke Universiteit Leuven
4 * Use of this software is governed by the MIT license
6 * Written by Sven Verdoolaege, K.U.Leuven, Departement
7 * Computerwetenschappen, Celestijnenlaan 200A, B-3001 Leuven, Belgium
8 */
10 #include <assert.h>
11 #include <string.h>
12 #include <strings.h>
13 #include <isl_map_private.h>
14 #include <isl/aff.h>
15 #include <isl/set.h>
16 #include "isl_tab.h"
17 #include "isl_sample.h"
18 #include "isl_scan.h"
19 #include <isl/seq.h>
20 #include <isl/ilp.h>
21 #include <isl/printer.h>
22 #include <isl_point_private.h>
23 #include <isl/options.h>
25 /* The input of this program is the same as that of the "example" program
26 * from the PipLib distribution, except that the "big parameter column"
27 * should always be -1.
29 * Context constraints in PolyLib format
30 * -1
31 * Problem constraints in PolyLib format
32 * Optional list of options
34 * The options are
35 * Maximize compute maximum instead of minimum
36 * Rational compute rational optimum instead of integer optimum
37 * Urs_parms don't assume parameters are non-negative
38 * Urs_unknowns don't assume unknowns are non-negative
41 struct options {
42 struct isl_options *isl;
43 unsigned verify;
44 unsigned format;
47 #define FORMAT_SET 0
48 #define FORMAT_AFF 1
50 struct isl_arg_choice pip_format[] = {
51 {"set", FORMAT_SET},
52 {"affine", FORMAT_AFF},
53 {0}
56 ISL_ARGS_START(struct options, options_args)
57 ISL_ARG_CHILD(struct options, isl, "isl", &isl_options_args, "isl options")
58 ISL_ARG_BOOL(struct options, verify, 'T', "verify", 0, NULL)
59 ISL_ARG_CHOICE(struct options, format, 0, "format",
60 pip_format, FORMAT_SET, "output format")
61 ISL_ARGS_END
63 ISL_ARG_DEF(options, struct options, options_args)
65 static __isl_give isl_basic_set *set_bounds(__isl_take isl_basic_set *bset)
67 unsigned nparam;
68 int i, r;
69 isl_point *pt, *pt2;
70 isl_basic_set *box;
72 nparam = isl_basic_set_dim(bset, isl_dim_param);
73 r = nparam >= 8 ? 4 : nparam >= 5 ? 6 : 30;
75 pt = isl_basic_set_sample_point(isl_basic_set_copy(bset));
76 pt2 = isl_point_copy(pt);
78 for (i = 0; i < nparam; ++i) {
79 pt = isl_point_add_ui(pt, isl_dim_param, i, r);
80 pt2 = isl_point_sub_ui(pt2, isl_dim_param, i, r);
83 box = isl_basic_set_box_from_points(pt, pt2);
85 return isl_basic_set_intersect(bset, box);
88 static struct isl_basic_set *to_parameter_domain(struct isl_basic_set *context)
90 context = isl_basic_set_move_dims(context, isl_dim_param, 0,
91 isl_dim_set, 0, isl_basic_set_dim(context, isl_dim_set));
92 context = isl_basic_set_params(context);
93 return context;
96 isl_basic_set *plug_in_parameters(isl_basic_set *bset, struct isl_vec *params)
98 int i;
100 for (i = 0; i < params->size - 1; ++i)
101 bset = isl_basic_set_fix(bset,
102 isl_dim_param, i, params->el[1 + i]);
104 bset = isl_basic_set_remove_dims(bset,
105 isl_dim_param, 0, params->size - 1);
107 isl_vec_free(params);
109 return bset;
112 isl_set *set_plug_in_parameters(isl_set *set, struct isl_vec *params)
114 int i;
116 for (i = 0; i < params->size - 1; ++i)
117 set = isl_set_fix(set, isl_dim_param, i, params->el[1 + i]);
119 set = isl_set_remove_dims(set, isl_dim_param, 0, params->size - 1);
121 isl_vec_free(params);
123 return set;
126 /* Compute the lexicographically minimal (or maximal if max is set)
127 * element of bset for the given values of the parameters, by
128 * successively solving an ilp problem in each direction.
130 struct isl_vec *opt_at(struct isl_basic_set *bset,
131 struct isl_vec *params, int max)
133 unsigned dim;
134 struct isl_vec *opt;
135 struct isl_vec *obj;
136 int i;
138 dim = isl_basic_set_dim(bset, isl_dim_set);
140 bset = plug_in_parameters(bset, params);
142 if (isl_basic_set_plain_is_empty(bset)) {
143 opt = isl_vec_alloc(bset->ctx, 0);
144 isl_basic_set_free(bset);
145 return opt;
148 opt = isl_vec_alloc(bset->ctx, 1 + dim);
149 assert(opt);
151 obj = isl_vec_alloc(bset->ctx, 1 + dim);
152 assert(obj);
154 isl_int_set_si(opt->el[0], 1);
155 isl_int_set_si(obj->el[0], 0);
157 for (i = 0; i < dim; ++i) {
158 enum isl_lp_result res;
160 isl_seq_clr(obj->el + 1, dim);
161 isl_int_set_si(obj->el[1 + i], 1);
162 res = isl_basic_set_solve_ilp(bset, max, obj->el,
163 &opt->el[1 + i], NULL);
164 if (res == isl_lp_empty)
165 goto empty;
166 assert(res == isl_lp_ok);
167 bset = isl_basic_set_fix(bset, isl_dim_set, i, opt->el[1 + i]);
170 isl_basic_set_free(bset);
171 isl_vec_free(obj);
173 return opt;
174 empty:
175 isl_vec_free(opt);
176 opt = isl_vec_alloc(bset->ctx, 0);
177 isl_basic_set_free(bset);
178 isl_vec_free(obj);
180 return opt;
183 struct isl_scan_pip {
184 struct isl_scan_callback callback;
185 isl_basic_set *bset;
186 isl_set *sol;
187 isl_set *empty;
188 int stride;
189 int n;
190 int max;
193 /* Check if the "manually" computed optimum of bset at the "sample"
194 * values of the parameters agrees with the solution of pilp problem
195 * represented by the pair (sol, empty).
196 * In particular, if there is no solution for this value of the parameters,
197 * then it should be an element of the parameter domain "empty".
198 * Otherwise, the optimal solution, should be equal to the result of
199 * plugging in the value of the parameters in "sol".
201 static int scan_one(struct isl_scan_callback *callback,
202 __isl_take isl_vec *sample)
204 struct isl_scan_pip *sp = (struct isl_scan_pip *)callback;
205 struct isl_vec *opt;
207 sp->n--;
209 opt = opt_at(isl_basic_set_copy(sp->bset), isl_vec_copy(sample), sp->max);
210 assert(opt);
212 if (opt->size == 0) {
213 isl_point *sample_pnt;
214 sample_pnt = isl_point_alloc(isl_set_get_space(sp->empty), sample);
215 assert(isl_set_contains_point(sp->empty, sample_pnt));
216 isl_point_free(sample_pnt);
217 isl_vec_free(opt);
218 } else {
219 isl_set *sol;
220 isl_set *opt_set;
221 opt_set = isl_set_from_basic_set(isl_basic_set_from_vec(opt));
222 sol = set_plug_in_parameters(isl_set_copy(sp->sol), sample);
223 assert(isl_set_is_equal(opt_set, sol));
224 isl_set_free(sol);
225 isl_set_free(opt_set);
228 if (!(sp->n % sp->stride)) {
229 printf("o");
230 fflush(stdout);
233 return sp->n >= 1 ? 0 : -1;
236 static void check_solution(isl_basic_set *bset, isl_basic_set *context,
237 isl_set *sol, isl_set *empty, int max)
239 struct isl_scan_pip sp;
240 isl_int count, count_max;
241 int i, n;
242 int r;
244 context = set_bounds(context);
245 context = isl_basic_set_underlying_set(context);
247 isl_int_init(count);
248 isl_int_init(count_max);
250 isl_int_set_si(count_max, 2000);
251 r = isl_basic_set_count_upto(context, count_max, &count);
252 assert(r >= 0);
253 n = isl_int_get_si(count);
255 isl_int_clear(count_max);
256 isl_int_clear(count);
258 sp.callback.add = scan_one;
259 sp.bset = bset;
260 sp.sol = sol;
261 sp.empty = empty;
262 sp.n = n;
263 sp.stride = n > 70 ? 1 + (n + 1)/70 : 1;
264 sp.max = max;
266 for (i = 0; i < n; i += sp.stride)
267 printf(".");
268 printf("\r");
269 fflush(stdout);
271 isl_basic_set_scan(context, &sp.callback);
273 printf("\n");
275 isl_basic_set_free(bset);
278 int main(int argc, char **argv)
280 struct isl_ctx *ctx;
281 struct isl_basic_set *context, *bset, *copy, *context_copy;
282 struct isl_set *set = NULL;
283 struct isl_set *empty;
284 isl_pw_multi_aff *pma = NULL;
285 int neg_one;
286 char s[1024];
287 int urs_parms = 0;
288 int urs_unknowns = 0;
289 int max = 0;
290 int rational = 0;
291 int n;
292 int nparam;
293 struct options *options;
295 options = options_new_with_defaults();
296 assert(options);
297 argc = options_parse(options, argc, argv, ISL_ARG_ALL);
299 ctx = isl_ctx_alloc_with_options(&options_args, options);
301 context = isl_basic_set_read_from_file(ctx, stdin);
302 assert(context);
303 n = fscanf(stdin, "%d", &neg_one);
304 assert(n == 1);
305 assert(neg_one == -1);
306 bset = isl_basic_set_read_from_file(ctx, stdin);
308 while (fgets(s, sizeof(s), stdin)) {
309 if (strncasecmp(s, "Maximize", 8) == 0)
310 max = 1;
311 if (strncasecmp(s, "Rational", 8) == 0) {
312 rational = 1;
313 bset = isl_basic_set_set_rational(bset);
315 if (strncasecmp(s, "Urs_parms", 9) == 0)
316 urs_parms = 1;
317 if (strncasecmp(s, "Urs_unknowns", 12) == 0)
318 urs_unknowns = 1;
320 if (!urs_parms)
321 context = isl_basic_set_intersect(context,
322 isl_basic_set_positive_orthant(isl_basic_set_get_space(context)));
323 context = to_parameter_domain(context);
324 nparam = isl_basic_set_dim(context, isl_dim_param);
325 if (nparam != isl_basic_set_dim(bset, isl_dim_param)) {
326 int dim = isl_basic_set_dim(bset, isl_dim_set);
327 bset = isl_basic_set_move_dims(bset, isl_dim_param, 0,
328 isl_dim_set, dim - nparam, nparam);
330 if (!urs_unknowns)
331 bset = isl_basic_set_intersect(bset,
332 isl_basic_set_positive_orthant(isl_basic_set_get_space(bset)));
334 if (options->verify) {
335 copy = isl_basic_set_copy(bset);
336 context_copy = isl_basic_set_copy(context);
339 if (options->format == FORMAT_AFF) {
340 if (max)
341 pma = isl_basic_set_partial_lexmax_pw_multi_aff(bset,
342 context, &empty);
343 else
344 pma = isl_basic_set_partial_lexmin_pw_multi_aff(bset,
345 context, &empty);
346 } else {
347 if (max)
348 set = isl_basic_set_partial_lexmax(bset,
349 context, &empty);
350 else
351 set = isl_basic_set_partial_lexmin(bset,
352 context, &empty);
355 if (options->verify) {
356 assert(!rational);
357 if (options->format == FORMAT_AFF)
358 set = isl_set_from_pw_multi_aff(pma);
359 check_solution(copy, context_copy, set, empty, max);
360 isl_set_free(set);
361 } else {
362 isl_printer *p;
363 p = isl_printer_to_file(ctx, stdout);
364 if (options->format == FORMAT_AFF)
365 p = isl_printer_print_pw_multi_aff(p, pma);
366 else
367 p = isl_printer_print_set(p, set);
368 p = isl_printer_end_line(p);
369 p = isl_printer_print_str(p, "no solution: ");
370 p = isl_printer_print_set(p, empty);
371 p = isl_printer_end_line(p);
372 isl_printer_free(p);
373 isl_set_free(set);
374 isl_pw_multi_aff_free(pma);
377 isl_set_free(empty);
378 isl_ctx_free(ctx);
380 return 0;