Line data Source code
1 : /*
2 : * Copyright 2010 INRIA Saclay
3 : *
4 : * Use of this software is governed by the MIT license
5 : *
6 : * Written by Sven Verdoolaege, INRIA Saclay - Ile-de-France,
7 : * Parc Club Orsay Universite, ZAC des vignes, 4 rue Jacques Monod,
8 : * 91893 Orsay, France
9 : */
10 :
11 : #include <isl_map_private.h>
12 : #include <isl/set.h>
13 : #include <isl_space_private.h>
14 : #include <isl_seq.h>
15 :
16 : /*
17 : * Let C be a cone and define
18 : *
19 : * C' := { y | forall x in C : y x >= 0 }
20 : *
21 : * C' contains the coefficients of all linear constraints
22 : * that are valid for C.
23 : * Furthermore, C'' = C.
24 : *
25 : * If C is defined as { x | A x >= 0 }
26 : * then any element in C' must be a non-negative combination
27 : * of the rows of A, i.e., y = t A with t >= 0. That is,
28 : *
29 : * C' = { y | exists t >= 0 : y = t A }
30 : *
31 : * If any of the rows in A actually represents an equality, then
32 : * also negative combinations of this row are allowed and so the
33 : * non-negativity constraint on the corresponding element of t
34 : * can be dropped.
35 : *
36 : * A polyhedron P = { x | b + A x >= 0 } can be represented
37 : * in homogeneous coordinates by the cone
38 : * C = { [z,x] | b z + A x >= and z >= 0 }
39 : * The valid linear constraints on C correspond to the valid affine
40 : * constraints on P.
41 : * This is essentially Farkas' lemma.
42 : *
43 : * Since
44 : * [ 1 0 ]
45 : * [ w y ] = [t_0 t] [ b A ]
46 : *
47 : * we have
48 : *
49 : * C' = { w, y | exists t_0, t >= 0 : y = t A and w = t_0 + t b }
50 : * or
51 : *
52 : * C' = { w, y | exists t >= 0 : y = t A and w - t b >= 0 }
53 : *
54 : * In practice, we introduce an extra variable (w), shifting all
55 : * other variables to the right, and an extra inequality
56 : * (w - t b >= 0) corresponding to the positivity constraint on
57 : * the homogeneous coordinate.
58 : *
59 : * When going back from coefficients to solutions, we immediately
60 : * plug in 1 for z, which corresponds to shifting all variables
61 : * to the left, with the leftmost ending up in the constant position.
62 : */
63 :
64 : /* Add the given prefix to all named isl_dim_set dimensions in "dim".
65 : */
66 0 : static __isl_give isl_space *isl_space_prefix(__isl_take isl_space *dim,
67 : const char *prefix)
68 : {
69 : int i;
70 : isl_ctx *ctx;
71 : unsigned nvar;
72 0 : size_t prefix_len = strlen(prefix);
73 :
74 0 : if (!dim)
75 0 : return NULL;
76 :
77 0 : ctx = isl_space_get_ctx(dim);
78 0 : nvar = isl_space_dim(dim, isl_dim_set);
79 :
80 0 : for (i = 0; i < nvar; ++i) {
81 : const char *name;
82 : char *prefix_name;
83 :
84 0 : name = isl_space_get_dim_name(dim, isl_dim_set, i);
85 0 : if (!name)
86 0 : continue;
87 :
88 0 : prefix_name = isl_alloc_array(ctx, char,
89 : prefix_len + strlen(name) + 1);
90 0 : if (!prefix_name)
91 0 : goto error;
92 0 : memcpy(prefix_name, prefix, prefix_len);
93 0 : strcpy(prefix_name + prefix_len, name);
94 :
95 0 : dim = isl_space_set_dim_name(dim, isl_dim_set, i, prefix_name);
96 0 : free(prefix_name);
97 : }
98 :
99 0 : return dim;
100 : error:
101 0 : isl_space_free(dim);
102 0 : return NULL;
103 : }
104 :
105 : /* Given a dimension specification of the solutions space, construct
106 : * a dimension specification for the space of coefficients.
107 : *
108 : * In particular transform
109 : *
110 : * [params] -> { S }
111 : *
112 : * to
113 : *
114 : * { coefficients[[cst, params] -> S] }
115 : *
116 : * and prefix each dimension name with "c_".
117 : */
118 0 : static __isl_give isl_space *isl_space_coefficients(__isl_take isl_space *dim)
119 : {
120 : isl_space *dim_param;
121 : unsigned nvar;
122 : unsigned nparam;
123 :
124 0 : nvar = isl_space_dim(dim, isl_dim_set);
125 0 : nparam = isl_space_dim(dim, isl_dim_param);
126 0 : dim_param = isl_space_copy(dim);
127 0 : dim_param = isl_space_drop_dims(dim_param, isl_dim_set, 0, nvar);
128 0 : dim_param = isl_space_move_dims(dim_param, isl_dim_set, 0,
129 : isl_dim_param, 0, nparam);
130 0 : dim_param = isl_space_prefix(dim_param, "c_");
131 0 : dim_param = isl_space_insert_dims(dim_param, isl_dim_set, 0, 1);
132 0 : dim_param = isl_space_set_dim_name(dim_param, isl_dim_set, 0, "c_cst");
133 0 : dim = isl_space_drop_dims(dim, isl_dim_param, 0, nparam);
134 0 : dim = isl_space_prefix(dim, "c_");
135 0 : dim = isl_space_join(isl_space_from_domain(dim_param),
136 : isl_space_from_range(dim));
137 0 : dim = isl_space_wrap(dim);
138 0 : dim = isl_space_set_tuple_name(dim, isl_dim_set, "coefficients");
139 :
140 0 : return dim;
141 : }
142 :
143 : /* Drop the given prefix from all named dimensions of type "type" in "dim".
144 : */
145 0 : static __isl_give isl_space *isl_space_unprefix(__isl_take isl_space *dim,
146 : enum isl_dim_type type, const char *prefix)
147 : {
148 : int i;
149 : unsigned n;
150 0 : size_t prefix_len = strlen(prefix);
151 :
152 0 : n = isl_space_dim(dim, type);
153 :
154 0 : for (i = 0; i < n; ++i) {
155 : const char *name;
156 :
157 0 : name = isl_space_get_dim_name(dim, type, i);
158 0 : if (!name)
159 0 : continue;
160 0 : if (strncmp(name, prefix, prefix_len))
161 0 : continue;
162 :
163 0 : dim = isl_space_set_dim_name(dim, type, i, name + prefix_len);
164 : }
165 :
166 0 : return dim;
167 : }
168 :
169 : /* Given a dimension specification of the space of coefficients, construct
170 : * a dimension specification for the space of solutions.
171 : *
172 : * In particular transform
173 : *
174 : * { coefficients[[cst, params] -> S] }
175 : *
176 : * to
177 : *
178 : * [params] -> { S }
179 : *
180 : * and drop the "c_" prefix from the dimension names.
181 : */
182 0 : static __isl_give isl_space *isl_space_solutions(__isl_take isl_space *dim)
183 : {
184 : unsigned nparam;
185 :
186 0 : dim = isl_space_unwrap(dim);
187 0 : dim = isl_space_drop_dims(dim, isl_dim_in, 0, 1);
188 0 : dim = isl_space_unprefix(dim, isl_dim_in, "c_");
189 0 : dim = isl_space_unprefix(dim, isl_dim_out, "c_");
190 0 : nparam = isl_space_dim(dim, isl_dim_in);
191 0 : dim = isl_space_move_dims(dim, isl_dim_param, 0, isl_dim_in, 0, nparam);
192 0 : dim = isl_space_range(dim);
193 :
194 0 : return dim;
195 : }
196 :
197 : /* Return the rational universe basic set in the given space.
198 : */
199 0 : static __isl_give isl_basic_set *rational_universe(__isl_take isl_space *space)
200 : {
201 : isl_basic_set *bset;
202 :
203 0 : bset = isl_basic_set_universe(space);
204 0 : bset = isl_basic_set_set_rational(bset);
205 :
206 0 : return bset;
207 : }
208 :
209 : /* Compute the dual of "bset" by applying Farkas' lemma.
210 : * As explained above, we add an extra dimension to represent
211 : * the coefficient of the constant term when going from solutions
212 : * to coefficients (shift == 1) and we drop the extra dimension when going
213 : * in the opposite direction (shift == -1). "dim" is the space in which
214 : * the dual should be created.
215 : *
216 : * If "bset" is (obviously) empty, then the way this emptiness
217 : * is represented by the constraints does not allow for the application
218 : * of the standard farkas algorithm. We therefore handle this case
219 : * specifically and return the universe basic set.
220 : */
221 0 : static __isl_give isl_basic_set *farkas(__isl_take isl_space *space,
222 : __isl_take isl_basic_set *bset, int shift)
223 : {
224 : int i, j, k;
225 0 : isl_basic_set *dual = NULL;
226 : unsigned total;
227 :
228 0 : if (isl_basic_set_plain_is_empty(bset)) {
229 0 : isl_basic_set_free(bset);
230 0 : return rational_universe(space);
231 : }
232 :
233 0 : total = isl_basic_set_total_dim(bset);
234 :
235 0 : dual = isl_basic_set_alloc_space(space, bset->n_eq + bset->n_ineq,
236 0 : total, bset->n_ineq + (shift > 0));
237 0 : dual = isl_basic_set_set_rational(dual);
238 :
239 0 : for (i = 0; i < bset->n_eq + bset->n_ineq; ++i) {
240 0 : k = isl_basic_set_alloc_div(dual);
241 0 : if (k < 0)
242 0 : goto error;
243 0 : isl_int_set_si(dual->div[k][0], 0);
244 : }
245 :
246 0 : for (i = 0; i < total; ++i) {
247 0 : k = isl_basic_set_alloc_equality(dual);
248 0 : if (k < 0)
249 0 : goto error;
250 0 : isl_seq_clr(dual->eq[k], 1 + shift + total);
251 0 : isl_int_set_si(dual->eq[k][1 + shift + i], -1);
252 0 : for (j = 0; j < bset->n_eq; ++j)
253 0 : isl_int_set(dual->eq[k][1 + shift + total + j],
254 : bset->eq[j][1 + i]);
255 0 : for (j = 0; j < bset->n_ineq; ++j)
256 0 : isl_int_set(dual->eq[k][1 + shift + total + bset->n_eq + j],
257 : bset->ineq[j][1 + i]);
258 : }
259 :
260 0 : for (i = 0; i < bset->n_ineq; ++i) {
261 0 : k = isl_basic_set_alloc_inequality(dual);
262 0 : if (k < 0)
263 0 : goto error;
264 0 : isl_seq_clr(dual->ineq[k],
265 0 : 1 + shift + total + bset->n_eq + bset->n_ineq);
266 0 : isl_int_set_si(dual->ineq[k][1 + shift + total + bset->n_eq + i], 1);
267 : }
268 :
269 0 : if (shift > 0) {
270 0 : k = isl_basic_set_alloc_inequality(dual);
271 0 : if (k < 0)
272 0 : goto error;
273 0 : isl_seq_clr(dual->ineq[k], 2 + total);
274 0 : isl_int_set_si(dual->ineq[k][1], 1);
275 0 : for (j = 0; j < bset->n_eq; ++j)
276 0 : isl_int_neg(dual->ineq[k][2 + total + j],
277 : bset->eq[j][0]);
278 0 : for (j = 0; j < bset->n_ineq; ++j)
279 0 : isl_int_neg(dual->ineq[k][2 + total + bset->n_eq + j],
280 : bset->ineq[j][0]);
281 : }
282 :
283 0 : dual = isl_basic_set_remove_divs(dual);
284 0 : dual = isl_basic_set_simplify(dual);
285 0 : dual = isl_basic_set_finalize(dual);
286 :
287 0 : isl_basic_set_free(bset);
288 0 : return dual;
289 : error:
290 0 : isl_basic_set_free(bset);
291 0 : isl_basic_set_free(dual);
292 0 : return NULL;
293 : }
294 :
295 : /* Construct a basic set containing the tuples of coefficients of all
296 : * valid affine constraints on the given basic set.
297 : */
298 0 : __isl_give isl_basic_set *isl_basic_set_coefficients(
299 : __isl_take isl_basic_set *bset)
300 : {
301 : isl_space *dim;
302 :
303 0 : if (!bset)
304 0 : return NULL;
305 0 : if (bset->n_div)
306 0 : isl_die(bset->ctx, isl_error_invalid,
307 : "input set not allowed to have local variables",
308 : goto error);
309 :
310 0 : dim = isl_basic_set_get_space(bset);
311 0 : dim = isl_space_coefficients(dim);
312 :
313 0 : return farkas(dim, bset, 1);
314 : error:
315 0 : isl_basic_set_free(bset);
316 0 : return NULL;
317 : }
318 :
319 : /* Construct a basic set containing the elements that satisfy all
320 : * affine constraints whose coefficient tuples are
321 : * contained in the given basic set.
322 : */
323 0 : __isl_give isl_basic_set *isl_basic_set_solutions(
324 : __isl_take isl_basic_set *bset)
325 : {
326 : isl_space *dim;
327 :
328 0 : if (!bset)
329 0 : return NULL;
330 0 : if (bset->n_div)
331 0 : isl_die(bset->ctx, isl_error_invalid,
332 : "input set not allowed to have local variables",
333 : goto error);
334 :
335 0 : dim = isl_basic_set_get_space(bset);
336 0 : dim = isl_space_solutions(dim);
337 :
338 0 : return farkas(dim, bset, -1);
339 : error:
340 0 : isl_basic_set_free(bset);
341 0 : return NULL;
342 : }
343 :
344 : /* Construct a basic set containing the tuples of coefficients of all
345 : * valid affine constraints on the given set.
346 : */
347 0 : __isl_give isl_basic_set *isl_set_coefficients(__isl_take isl_set *set)
348 : {
349 : int i;
350 : isl_basic_set *coeff;
351 :
352 0 : if (!set)
353 0 : return NULL;
354 0 : if (set->n == 0) {
355 0 : isl_space *space = isl_set_get_space(set);
356 0 : space = isl_space_coefficients(space);
357 0 : isl_set_free(set);
358 0 : return rational_universe(space);
359 : }
360 :
361 0 : coeff = isl_basic_set_coefficients(isl_basic_set_copy(set->p[0]));
362 :
363 0 : for (i = 1; i < set->n; ++i) {
364 : isl_basic_set *bset, *coeff_i;
365 0 : bset = isl_basic_set_copy(set->p[i]);
366 0 : coeff_i = isl_basic_set_coefficients(bset);
367 0 : coeff = isl_basic_set_intersect(coeff, coeff_i);
368 : }
369 :
370 0 : isl_set_free(set);
371 0 : return coeff;
372 : }
373 :
374 : /* Wrapper around isl_basic_set_coefficients for use
375 : * as a isl_basic_set_list_map callback.
376 : */
377 0 : static __isl_give isl_basic_set *coefficients_wrap(
378 : __isl_take isl_basic_set *bset, void *user)
379 : {
380 0 : return isl_basic_set_coefficients(bset);
381 : }
382 :
383 : /* Replace the elements of "list" by the result of applying
384 : * isl_basic_set_coefficients to them.
385 : */
386 0 : __isl_give isl_basic_set_list *isl_basic_set_list_coefficients(
387 : __isl_take isl_basic_set_list *list)
388 : {
389 0 : return isl_basic_set_list_map(list, &coefficients_wrap, NULL);
390 : }
391 :
392 : /* Construct a basic set containing the elements that satisfy all
393 : * affine constraints whose coefficient tuples are
394 : * contained in the given set.
395 : */
396 0 : __isl_give isl_basic_set *isl_set_solutions(__isl_take isl_set *set)
397 : {
398 : int i;
399 : isl_basic_set *sol;
400 :
401 0 : if (!set)
402 0 : return NULL;
403 0 : if (set->n == 0) {
404 0 : isl_space *space = isl_set_get_space(set);
405 0 : space = isl_space_solutions(space);
406 0 : isl_set_free(set);
407 0 : return rational_universe(space);
408 : }
409 :
410 0 : sol = isl_basic_set_solutions(isl_basic_set_copy(set->p[0]));
411 :
412 0 : for (i = 1; i < set->n; ++i) {
413 : isl_basic_set *bset, *sol_i;
414 0 : bset = isl_basic_set_copy(set->p[i]);
415 0 : sol_i = isl_basic_set_solutions(bset);
416 0 : sol = isl_basic_set_intersect(sol, sol_i);
417 : }
418 :
419 0 : isl_set_free(set);
420 0 : return sol;
421 : }
|
