1/*
2 * Copyright 2008-2009 Katholieke Universiteit Leuven
3 *
4 * Use of this software is governed by the MIT 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/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>
24
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.
28 *
29 * Context constraints in PolyLib format
30 * -1
31 * Problem constraints in PolyLib format
32 * Optional list of options
33 *
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
39 */
40
41struct options {
42	struct isl_options	*isl;
43	unsigned		 verify;
44	unsigned		 format;
45};
46
47#define FORMAT_SET	0
48#define FORMAT_AFF	1
49
50struct isl_arg_choice pip_format[] = {
51	{"set",		FORMAT_SET},
52	{"affine",	FORMAT_AFF},
53	{0}
54};
55
56ISL_ARGS_START(struct options, options_args)
57ISL_ARG_CHILD(struct options, isl, "isl", &isl_options_args, "isl options")
58ISL_ARG_BOOL(struct options, verify, 'T', "verify", 0, NULL)
59ISL_ARG_CHOICE(struct options, format, 0, "format",
60	pip_format, FORMAT_SET, "output format")
61ISL_ARGS_END
62
63ISL_ARG_DEF(options, struct options, options_args)
64
65static __isl_give isl_basic_set *set_bounds(__isl_take isl_basic_set *bset)
66{
67	unsigned nparam;
68	int i, r;
69	isl_point *pt, *pt2;
70	isl_basic_set *box;
71
72	nparam = isl_basic_set_dim(bset, isl_dim_param);
73	r = nparam >= 8 ? 4 : nparam >= 5 ? 6 : 30;
74
75	pt = isl_basic_set_sample_point(isl_basic_set_copy(bset));
76	pt2 = isl_point_copy(pt);
77
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);
81	}
82
83	box = isl_basic_set_box_from_points(pt, pt2);
84
85	return isl_basic_set_intersect(bset, box);
86}
87
88static struct isl_basic_set *to_parameter_domain(struct isl_basic_set *context)
89{
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;
94}
95
96isl_basic_set *plug_in_parameters(isl_basic_set *bset, struct isl_vec *params)
97{
98	int i;
99
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]);
103
104	bset = isl_basic_set_remove_dims(bset,
105					 isl_dim_param, 0, params->size - 1);
106
107	isl_vec_free(params);
108
109	return bset;
110}
111
112isl_set *set_plug_in_parameters(isl_set *set, struct isl_vec *params)
113{
114	int i;
115
116	for (i = 0; i < params->size - 1; ++i)
117		set = isl_set_fix(set, isl_dim_param, i, params->el[1 + i]);
118
119	set = isl_set_remove_dims(set, isl_dim_param, 0, params->size - 1);
120
121	isl_vec_free(params);
122
123	return set;
124}
125
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.
129 */
130struct isl_vec *opt_at(struct isl_basic_set *bset,
131	struct isl_vec *params, int max)
132{
133	unsigned dim;
134	struct isl_vec *opt;
135	struct isl_vec *obj;
136	int i;
137
138	dim = isl_basic_set_dim(bset, isl_dim_set);
139
140	bset = plug_in_parameters(bset, params);
141
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;
146	}
147
148	opt = isl_vec_alloc(bset->ctx, 1 + dim);
149	assert(opt);
150
151	obj = isl_vec_alloc(bset->ctx, 1 + dim);
152	assert(obj);
153
154	isl_int_set_si(opt->el[0], 1);
155	isl_int_set_si(obj->el[0], 0);
156
157	for (i = 0; i < dim; ++i) {
158		enum isl_lp_result res;
159
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]);
168	}
169
170	isl_basic_set_free(bset);
171	isl_vec_free(obj);
172
173	return opt;
174empty:
175	isl_vec_free(opt);
176	opt = isl_vec_alloc(bset->ctx, 0);
177	isl_basic_set_free(bset);
178	isl_vec_free(obj);
179
180	return opt;
181}
182
183struct 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;
191};
192
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".
200 */
201static int scan_one(struct isl_scan_callback *callback,
202	__isl_take isl_vec *sample)
203{
204	struct isl_scan_pip *sp = (struct isl_scan_pip *)callback;
205	struct isl_vec *opt;
206
207	sp->n--;
208
209	opt = opt_at(isl_basic_set_copy(sp->bset), isl_vec_copy(sample), sp->max);
210	assert(opt);
211
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);
226	}
227
228	if (!(sp->n % sp->stride)) {
229		printf("o");
230		fflush(stdout);
231	}
232
233	return sp->n >= 1 ? 0 : -1;
234}
235
236static void check_solution(isl_basic_set *bset, isl_basic_set *context,
237	isl_set *sol, isl_set *empty, int max)
238{
239	struct isl_scan_pip sp;
240	isl_int count, count_max;
241	int i, n;
242	int r;
243
244	context = set_bounds(context);
245	context = isl_basic_set_underlying_set(context);
246
247	isl_int_init(count);
248	isl_int_init(count_max);
249
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);
254
255	isl_int_clear(count_max);
256	isl_int_clear(count);
257
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;
265
266	for (i = 0; i < n; i += sp.stride)
267		printf(".");
268	printf("\r");
269	fflush(stdout);
270
271	isl_basic_set_scan(context, &sp.callback);
272
273	printf("\n");
274
275	isl_basic_set_free(bset);
276}
277
278int main(int argc, char **argv)
279{
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;
294
295	options = options_new_with_defaults();
296	assert(options);
297	argc = options_parse(options, argc, argv, ISL_ARG_ALL);
298
299	ctx = isl_ctx_alloc_with_options(&options_args, options);
300
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);
307
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);
314		}
315		if (strncasecmp(s, "Urs_parms", 9) == 0)
316			urs_parms = 1;
317		if (strncasecmp(s, "Urs_unknowns", 12) == 0)
318			urs_unknowns = 1;
319	}
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);
329	}
330	if (!urs_unknowns)
331		bset = isl_basic_set_intersect(bset,
332		isl_basic_set_positive_orthant(isl_basic_set_get_space(bset)));
333
334	if (options->verify) {
335		copy = isl_basic_set_copy(bset);
336		context_copy = isl_basic_set_copy(context);
337	}
338
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);
353	}
354
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);
375	}
376
377	isl_set_free(empty);
378	isl_ctx_free(ctx);
379
380	return 0;
381}
382