add packaging
[platform/upstream/isl.git] / isl_map_piplib.c
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 <isl_ctx_private.h>
11 #include <isl_map_private.h>
12 #include <isl/set.h>
13 #include <isl/map.h>
14 #include <isl/mat.h>
15 #include <isl/seq.h>
16 #include "isl_piplib.h"
17 #include "isl_map_piplib.h"
18
19 static void copy_values_from(isl_int *dst, Entier *src, unsigned n)
20 {
21         int i;
22
23         for (i = 0; i < n; ++i)
24                 entier_assign(dst[i], src[i]);
25 }
26
27 static void add_value(isl_int *dst, Entier *src)
28 {
29         mpz_add(*dst, *dst, *src);
30 }
31
32 static void copy_constraint_from(isl_int *dst, PipVector *src,
33                 unsigned nparam, unsigned n_in, unsigned n_out,
34                 unsigned extra, int *pos)
35 {
36         int i;
37
38         copy_values_from(dst, src->the_vector+src->nb_elements-1, 1);
39         copy_values_from(dst+1, src->the_vector, nparam+n_in);
40         isl_seq_clr(dst+1+nparam+n_in, n_out);
41         isl_seq_clr(dst+1+nparam+n_in+n_out, extra);
42         for (i = 0; i + n_in + nparam < src->nb_elements-1; ++i) {
43                 int p = pos[i];
44                 add_value(&dst[1+nparam+n_in+n_out+p],
45                           &src->the_vector[n_in+nparam+i]);
46         }
47 }
48
49 static int add_inequality(struct isl_ctx *ctx,
50                    struct isl_basic_map *bmap, int *pos, PipVector *vec)
51 {
52         unsigned nparam = isl_basic_map_n_param(bmap);
53         unsigned n_in = isl_basic_map_n_in(bmap);
54         unsigned n_out = isl_basic_map_n_out(bmap);
55         unsigned n_div = isl_basic_map_n_div(bmap);
56         int i = isl_basic_map_alloc_inequality(bmap);
57         if (i < 0)
58                 return -1;
59         copy_constraint_from(bmap->ineq[i], vec,
60             nparam, n_in, n_out, n_div, pos);
61
62         return i;
63 }
64
65 /* For a div d = floor(f/m), add the constraints
66  *
67  *              f - m d >= 0
68  *              -(f-(n-1)) + m d >= 0
69  *
70  * Note that the second constraint is the negation of
71  *
72  *              f - m d >= n
73  */
74 static int add_div_constraints(struct isl_ctx *ctx,
75         struct isl_basic_map *bmap, int *pos, PipNewparm *p, unsigned div)
76 {
77         int i, j;
78         unsigned total = isl_basic_map_total_dim(bmap);
79         unsigned div_pos = 1 + total - bmap->n_div + div;
80
81         i = add_inequality(ctx, bmap, pos, p->vector);
82         if (i < 0)
83                 return -1;
84         copy_values_from(&bmap->ineq[i][div_pos], &p->deno, 1);
85         isl_int_neg(bmap->ineq[i][div_pos], bmap->ineq[i][div_pos]);
86
87         j = isl_basic_map_alloc_inequality(bmap);
88         if (j < 0)
89                 return -1;
90         isl_seq_neg(bmap->ineq[j], bmap->ineq[i], 1 + total);
91         isl_int_add(bmap->ineq[j][0], bmap->ineq[j][0], bmap->ineq[j][div_pos]);
92         isl_int_sub_ui(bmap->ineq[j][0], bmap->ineq[j][0], 1);
93         return j;
94 }
95
96 static int add_equality(struct isl_ctx *ctx,
97                    struct isl_basic_map *bmap, int *pos,
98                    unsigned var, PipVector *vec)
99 {
100         int i;
101         unsigned nparam = isl_basic_map_n_param(bmap);
102         unsigned n_in = isl_basic_map_n_in(bmap);
103         unsigned n_out = isl_basic_map_n_out(bmap);
104
105         isl_assert(ctx, var < n_out, return -1);
106
107         i = isl_basic_map_alloc_equality(bmap);
108         if (i < 0)
109                 return -1;
110         copy_constraint_from(bmap->eq[i], vec,
111             nparam, n_in, n_out, bmap->extra, pos);
112         isl_int_set_si(bmap->eq[i][1+nparam+n_in+var], -1);
113
114         return i;
115 }
116
117 static int find_div(struct isl_ctx *ctx,
118                    struct isl_basic_map *bmap, int *pos, PipNewparm *p)
119 {
120         int i, j;
121         unsigned nparam = isl_basic_map_n_param(bmap);
122         unsigned n_in = isl_basic_map_n_in(bmap);
123         unsigned n_out = isl_basic_map_n_out(bmap);
124
125         i = isl_basic_map_alloc_div(bmap);
126         if (i < 0)
127                 return -1;
128
129         copy_constraint_from(bmap->div[i]+1, p->vector,
130             nparam, n_in, n_out, bmap->extra, pos);
131
132         copy_values_from(bmap->div[i], &p->deno, 1);
133         for (j = 0; j < i; ++j)
134                 if (isl_seq_eq(bmap->div[i], bmap->div[j],
135                                 1+1+isl_basic_map_total_dim(bmap)+j)) {
136                         isl_basic_map_free_div(bmap, 1);
137                         return j;
138                 }
139
140         if (add_div_constraints(ctx, bmap, pos, p, i) < 0)
141                 return -1;
142
143         return i;
144 }
145
146 /* Count some properties of a quast
147  * - maximal number of new parameters
148  * - maximal depth
149  * - total number of solutions
150  * - total number of empty branches
151  */
152 static void quast_count(PipQuast *q, int *maxnew, int depth, int *maxdepth,
153                         int *sol, int *nosol)
154 {
155         PipNewparm *p;
156
157         for (p = q->newparm; p; p = p->next)
158                 if (p->rank > *maxnew)
159                         *maxnew = p->rank;
160         if (q->condition) {
161                 if (++depth > *maxdepth)
162                         *maxdepth = depth;
163                 quast_count(q->next_else, maxnew, depth, maxdepth, sol, nosol);
164                 quast_count(q->next_then, maxnew, depth, maxdepth, sol, nosol);
165         } else {
166                 if (q->list)
167                         ++(*sol);
168                 else
169                         ++(*nosol);
170         }
171 }
172
173 /*
174  * pos: array of length bmap->set.extra, mapping each of the existential
175  *              variables PIP proposes to an existential variable in bmap
176  * bmap: collects the currently active constraints
177  * rest: collects the empty leaves of the quast (if not NULL)
178  */
179 struct scan_data {
180         struct isl_ctx                  *ctx;
181         struct isl_basic_map            *bmap;
182         struct isl_set                  **rest;
183         int        *pos;
184 };
185
186 /*
187  * New existentially quantified variables are places after the existing ones.
188  */
189 static struct isl_map *scan_quast_r(struct scan_data *data, PipQuast *q,
190                                     struct isl_map *map)
191 {
192         PipNewparm *p;
193         struct isl_basic_map *bmap = data->bmap;
194         unsigned old_n_div = bmap->n_div;
195         unsigned nparam = isl_basic_map_n_param(bmap);
196         unsigned n_in = isl_basic_map_n_in(bmap);
197         unsigned n_out = isl_basic_map_n_out(bmap);
198
199         if (!map)
200                 goto error;
201
202         for (p = q->newparm; p; p = p->next) {
203                 int pos;
204                 unsigned pip_param = nparam + n_in;
205
206                 pos = find_div(data->ctx, bmap, data->pos, p);
207                 if (pos < 0)
208                         goto error;
209                 data->pos[p->rank - pip_param] = pos;
210         }
211
212         if (q->condition) {
213                 int pos = add_inequality(data->ctx, bmap, data->pos,
214                                          q->condition);
215                 if (pos < 0)
216                         goto error;
217                 map = scan_quast_r(data, q->next_then, map);
218
219                 if (isl_inequality_negate(bmap, pos))
220                         goto error;
221                 map = scan_quast_r(data, q->next_else, map);
222
223                 if (isl_basic_map_free_inequality(bmap, 1))
224                         goto error;
225         } else if (q->list) {
226                 PipList *l;
227                 int j;
228                 /* if bmap->n_out is zero, we are only interested in the domains
229                  * where a solution exists and not in the actual solution
230                  */
231                 for (j = 0, l = q->list; j < n_out && l; ++j, l = l->next)
232                         if (add_equality(data->ctx, bmap, data->pos, j,
233                                                 l->vector) < 0)
234                                 goto error;
235                 map = isl_map_add_basic_map(map, isl_basic_map_copy(bmap));
236                 if (isl_basic_map_free_equality(bmap, n_out))
237                         goto error;
238         } else if (data->rest) {
239                 struct isl_basic_set *bset;
240                 bset = isl_basic_set_from_basic_map(isl_basic_map_copy(bmap));
241                 bset = isl_basic_set_drop_dims(bset, n_in, n_out);
242                 if (!bset)
243                         goto error;
244                 *data->rest = isl_set_add_basic_set(*data->rest, bset);
245         }
246
247         if (isl_basic_map_free_inequality(bmap, 2*(bmap->n_div - old_n_div)))
248                 goto error;
249         if (isl_basic_map_free_div(bmap, bmap->n_div - old_n_div))
250                 goto error;
251         return map;
252 error:
253         isl_map_free(map);
254         return NULL;
255 }
256
257 /*
258  * Returns a map of dimension "keep_dim" with "context" as domain and
259  * as range the first "isl_space_dim(keep_dim, isl_dim_out)" variables
260  * in the quast lists.
261  */
262 static struct isl_map *isl_map_from_quast(struct isl_ctx *ctx, PipQuast *q,
263                 isl_space *keep_dim,
264                 struct isl_basic_set *context,
265                 struct isl_set **rest)
266 {
267         int             pip_param;
268         int             nexist;
269         int             max_depth;
270         int             n_sol, n_nosol;
271         struct scan_data        data;
272         struct isl_map          *map = NULL;
273         isl_space               *dims;
274         unsigned                nparam;
275         unsigned                dim;
276         unsigned                keep;
277
278         data.ctx = ctx;
279         data.rest = rest;
280         data.bmap = NULL;
281         data.pos = NULL;
282
283         if (!context || !keep_dim)
284                 goto error;
285
286         dim = isl_basic_set_n_dim(context);
287         nparam = isl_basic_set_n_param(context);
288         keep = isl_space_dim(keep_dim, isl_dim_out);
289         pip_param = nparam + dim;
290
291         max_depth = 0;
292         n_sol = 0;
293         n_nosol = 0;
294         nexist = pip_param-1;
295         quast_count(q, &nexist, 0, &max_depth, &n_sol, &n_nosol);
296         nexist -= pip_param-1;
297
298         if (rest) {
299                 *rest = isl_set_alloc_space(isl_space_copy(context->dim), n_nosol,
300                                         ISL_MAP_DISJOINT);
301                 if (!*rest)
302                         goto error;
303         }
304         map = isl_map_alloc_space(isl_space_copy(keep_dim), n_sol,
305                                 ISL_MAP_DISJOINT);
306         if (!map)
307                 goto error;
308
309         dims = isl_space_reverse(isl_space_copy(context->dim));
310         data.bmap = isl_basic_map_from_basic_set(context, dims);
311         data.bmap = isl_basic_map_extend_space(data.bmap,
312                 keep_dim, nexist, keep, max_depth+2*nexist);
313         if (!data.bmap)
314                 goto error2;
315
316         if (data.bmap->extra) {
317                 int i;
318                 data.pos = isl_alloc_array(ctx, int, data.bmap->extra);
319                 if (!data.pos)
320                         goto error;
321                 for (i = 0; i < data.bmap->n_div; ++i)
322                         data.pos[i] = i;
323         }
324
325         map = scan_quast_r(&data, q, map);
326         map = isl_map_finalize(map);
327         if (!map)
328                 goto error2;
329         if (rest) {
330                 *rest = isl_set_finalize(*rest);
331                 if (!*rest)
332                         goto error2;
333         }
334         isl_basic_map_free(data.bmap);
335         if (data.pos)
336                 free(data.pos);
337         return map;
338 error:
339         isl_basic_set_free(context);
340         isl_space_free(keep_dim);
341 error2:
342         if (data.pos)
343                 free(data.pos);
344         isl_basic_map_free(data.bmap);
345         isl_map_free(map);
346         if (rest) {
347                 isl_set_free(*rest);
348                 *rest = NULL;
349         }
350         return NULL;
351 }
352
353 static void copy_values_to(Entier *dst, isl_int *src, unsigned n)
354 {
355         int i;
356
357         for (i = 0; i < n; ++i)
358                 entier_assign(dst[i], src[i]);
359 }
360
361 static void copy_constraint_to(Entier *dst, isl_int *src,
362                 unsigned pip_param, unsigned pip_var,
363                 unsigned extra_front, unsigned extra_back)
364 {
365         copy_values_to(dst+1+extra_front+pip_var+pip_param+extra_back, src, 1);
366         copy_values_to(dst+1+extra_front+pip_var, src+1, pip_param);
367         copy_values_to(dst+1+extra_front, src+1+pip_param, pip_var);
368 }
369
370 PipMatrix *isl_basic_map_to_pip(struct isl_basic_map *bmap, unsigned pip_param,
371                          unsigned extra_front, unsigned extra_back)
372 {
373         int i;
374         unsigned nrow;
375         unsigned ncol;
376         PipMatrix *M;
377         unsigned off;
378         unsigned pip_var = isl_basic_map_total_dim(bmap) - pip_param;
379
380         nrow = extra_front + bmap->n_eq + bmap->n_ineq;
381         ncol = 1 + extra_front + pip_var + pip_param + extra_back + 1;
382         M = pip_matrix_alloc(nrow, ncol);
383         if (!M)
384                 return NULL;
385
386         off = extra_front;
387         for (i = 0; i < bmap->n_eq; ++i) {
388                 entier_set_si(M->p[off+i][0], 0);
389                 copy_constraint_to(M->p[off+i], bmap->eq[i],
390                                    pip_param, pip_var, extra_front, extra_back);
391         }
392         off += bmap->n_eq;
393         for (i = 0; i < bmap->n_ineq; ++i) {
394                 entier_set_si(M->p[off+i][0], 1);
395                 copy_constraint_to(M->p[off+i], bmap->ineq[i],
396                                    pip_param, pip_var, extra_front, extra_back);
397         }
398         return M;
399 }
400
401 PipMatrix *isl_basic_set_to_pip(struct isl_basic_set *bset, unsigned pip_param,
402                          unsigned extra_front, unsigned extra_back)
403 {
404         return isl_basic_map_to_pip((struct isl_basic_map *)bset,
405                                         pip_param, extra_front, extra_back);
406 }
407
408 struct isl_map *isl_pip_basic_map_lexopt(
409                 struct isl_basic_map *bmap, struct isl_basic_set *dom,
410                 struct isl_set **empty, int max)
411 {
412         PipOptions      *options;
413         PipQuast        *sol;
414         struct isl_map  *map;
415         struct isl_ctx  *ctx;
416         PipMatrix *domain = NULL, *context = NULL;
417         unsigned         nparam, n_in, n_out;
418
419         bmap = isl_basic_map_detect_equalities(bmap);
420         if (!bmap || !dom)
421                 goto error;
422
423         ctx = bmap->ctx;
424         isl_assert(ctx, isl_basic_map_compatible_domain(bmap, dom), goto error);
425         nparam = isl_basic_map_n_param(bmap);
426         n_in = isl_basic_map_n_in(bmap);
427         n_out = isl_basic_map_n_out(bmap);
428
429         domain = isl_basic_map_to_pip(bmap, nparam + n_in, 0, dom->n_div);
430         if (!domain)
431                 goto error;
432         context = isl_basic_map_to_pip((struct isl_basic_map *)dom, 0, 0, 0);
433         if (!context)
434                 goto error;
435
436         options = pip_options_init();
437         options->Simplify = 1;
438         options->Maximize = max;
439         options->Urs_unknowns = -1;
440         options->Urs_parms = -1;
441         sol = pip_solve(domain, context, -1, options);
442
443         if (sol) {
444                 struct isl_basic_set *copy;
445                 copy = isl_basic_set_copy(dom);
446                 map = isl_map_from_quast(ctx, sol,
447                                 isl_space_copy(bmap->dim), copy, empty);
448         } else {
449                 map = isl_map_empty_like_basic_map(bmap);
450                 if (empty)
451                         *empty = NULL;
452         }
453         if (!map)
454                 goto error;
455         if (map->n == 0 && empty) {
456                 isl_set_free(*empty);
457                 *empty = isl_set_from_basic_set(dom);
458         } else
459                 isl_basic_set_free(dom);
460         isl_basic_map_free(bmap);
461
462         pip_quast_free(sol);
463         pip_options_free(options);
464         pip_matrix_free(domain);
465         pip_matrix_free(context);
466
467         return map;
468 error:
469         if (domain)
470                 pip_matrix_free(domain);
471         if (context)
472                 pip_matrix_free(context);
473         isl_basic_map_free(bmap);
474         isl_basic_set_free(dom);
475         return NULL;
476 }