+isl_basic_set *plug_in_parameters(isl_basic_set *bset, struct isl_vec *params)
+{
+ int i;
+
+ for (i = 0; i < params->size - 1; ++i)
+ bset = isl_basic_set_fix(bset,
+ isl_dim_param, i, params->el[1 + i]);
+
+ bset = isl_basic_set_remove_dims(bset,
+ isl_dim_param, 0, params->size - 1);
+
+ isl_vec_free(params);
+
+ return bset;
+}
+
+isl_set *set_plug_in_parameters(isl_set *set, struct isl_vec *params)
+{
+ int i;
+
+ for (i = 0; i < params->size - 1; ++i)
+ set = isl_set_fix(set, isl_dim_param, i, params->el[1 + i]);
+
+ set = isl_set_remove_dims(set, isl_dim_param, 0, params->size - 1);
+
+ isl_vec_free(params);
+
+ return set;
+}
+
+/* Compute the lexicographically minimal (or maximal if max is set)
+ * element of bset for the given values of the parameters, by
+ * successively solving an ilp problem in each direction.
+ */
+struct isl_vec *opt_at(struct isl_basic_set *bset,
+ struct isl_vec *params, int max)
+{
+ unsigned dim;
+ struct isl_vec *opt;
+ struct isl_vec *obj;
+ int i;
+
+ dim = isl_basic_set_dim(bset, isl_dim_set);
+
+ bset = plug_in_parameters(bset, params);
+
+ if (isl_basic_set_plain_is_empty(bset)) {
+ opt = isl_vec_alloc(bset->ctx, 0);
+ isl_basic_set_free(bset);
+ return opt;
+ }
+
+ opt = isl_vec_alloc(bset->ctx, 1 + dim);
+ assert(opt);
+
+ obj = isl_vec_alloc(bset->ctx, 1 + dim);
+ assert(obj);
+
+ isl_int_set_si(opt->el[0], 1);
+ isl_int_set_si(obj->el[0], 0);
+
+ for (i = 0; i < dim; ++i) {
+ enum isl_lp_result res;
+
+ isl_seq_clr(obj->el + 1, dim);
+ isl_int_set_si(obj->el[1 + i], 1);
+ res = isl_basic_set_solve_ilp(bset, max, obj->el,
+ &opt->el[1 + i], NULL);
+ if (res == isl_lp_empty)
+ goto empty;
+ assert(res == isl_lp_ok);
+ bset = isl_basic_set_fix(bset, isl_dim_set, i, opt->el[1 + i]);
+ }
+
+ isl_basic_set_free(bset);
+ isl_vec_free(obj);
+
+ return opt;
+empty:
+ isl_vec_free(opt);
+ opt = isl_vec_alloc(bset->ctx, 0);
+ isl_basic_set_free(bset);
+ isl_vec_free(obj);
+
+ return opt;
+}
+
+struct isl_scan_pip {
+ struct isl_scan_callback callback;
+ isl_basic_set *bset;
+ isl_set *sol;
+ isl_set *empty;
+ int stride;
+ int n;
+ int max;
+};
+
+/* Check if the "manually" computed optimum of bset at the "sample"
+ * values of the parameters agrees with the solution of pilp problem
+ * represented by the pair (sol, empty).
+ * In particular, if there is no solution for this value of the parameters,
+ * then it should be an element of the parameter domain "empty".
+ * Otherwise, the optimal solution, should be equal to the result of
+ * plugging in the value of the parameters in "sol".
+ */
+static int scan_one(struct isl_scan_callback *callback,
+ __isl_take isl_vec *sample)
+{
+ struct isl_scan_pip *sp = (struct isl_scan_pip *)callback;
+ struct isl_vec *opt;
+
+ sp->n--;
+
+ opt = opt_at(isl_basic_set_copy(sp->bset), isl_vec_copy(sample), sp->max);
+ assert(opt);
+
+ if (opt->size == 0) {
+ isl_point *sample_pnt;
+ sample_pnt = isl_point_alloc(isl_set_get_space(sp->empty), sample);
+ assert(isl_set_contains_point(sp->empty, sample_pnt));
+ isl_point_free(sample_pnt);
+ isl_vec_free(opt);
+ } else {
+ isl_set *sol;
+ isl_set *opt_set;
+ opt_set = isl_set_from_basic_set(isl_basic_set_from_vec(opt));
+ sol = set_plug_in_parameters(isl_set_copy(sp->sol), sample);
+ assert(isl_set_is_equal(opt_set, sol));
+ isl_set_free(sol);
+ isl_set_free(opt_set);
+ }
+
+ if (!(sp->n % sp->stride)) {
+ printf("o");
+ fflush(stdout);
+ }
+
+ return sp->n >= 1 ? 0 : -1;
+}
+
+static void check_solution(isl_basic_set *bset, isl_basic_set *context,
+ isl_set *sol, isl_set *empty, int max)
+{
+ struct isl_scan_pip sp;
+ isl_int count, count_max;
+ int i, n;
+ int r;
+
+ context = set_bounds(context);
+ context = isl_basic_set_underlying_set(context);
+
+ isl_int_init(count);
+ isl_int_init(count_max);
+
+ isl_int_set_si(count_max, 2000);
+ r = isl_basic_set_count_upto(context, count_max, &count);
+ assert(r >= 0);
+ n = isl_int_get_si(count);
+
+ isl_int_clear(count_max);
+ isl_int_clear(count);
+
+ sp.callback.add = scan_one;
+ sp.bset = bset;
+ sp.sol = sol;
+ sp.empty = empty;
+ sp.n = n;
+ sp.stride = n > 70 ? 1 + (n + 1)/70 : 1;
+ sp.max = max;
+
+ for (i = 0; i < n; i += sp.stride)
+ printf(".");
+ printf("\r");
+ fflush(stdout);
+
+ isl_basic_set_scan(context, &sp.callback);
+
+ printf("\n");
+
+ isl_basic_set_free(bset);
+}
+