add isl_aff_add_constant_num{,_si}
[platform/upstream/isl.git] / isl_aff.c
1 /*
2  * Copyright 2011      INRIA Saclay
3  * Copyright 2011      Sven Verdoolaege
4  * Copyright 2012      Ecole Normale Superieure
5  *
6  * Use of this software is governed by the GNU LGPLv2.1 license
7  *
8  * Written by Sven Verdoolaege, INRIA Saclay - Ile-de-France,
9  * Parc Club Orsay Universite, ZAC des vignes, 4 rue Jacques Monod,
10  * 91893 Orsay, France
11  * and Ecole Normale Superieure, 45 rue d’Ulm, 75230 Paris, France
12  */
13
14 #include <isl_ctx_private.h>
15 #define ISL_DIM_H
16 #include <isl_map_private.h>
17 #include <isl_union_map_private.h>
18 #include <isl_aff_private.h>
19 #include <isl_space_private.h>
20 #include <isl_local_space_private.h>
21 #include <isl_mat_private.h>
22 #include <isl_list_private.h>
23 #include <isl/constraint.h>
24 #include <isl/seq.h>
25 #include <isl/set.h>
26 #include <isl_config.h>
27
28 __isl_give isl_aff *isl_aff_alloc_vec(__isl_take isl_local_space *ls,
29         __isl_take isl_vec *v)
30 {
31         isl_aff *aff;
32
33         if (!ls || !v)
34                 goto error;
35
36         aff = isl_calloc_type(v->ctx, struct isl_aff);
37         if (!aff)
38                 goto error;
39
40         aff->ref = 1;
41         aff->ls = ls;
42         aff->v = v;
43
44         return aff;
45 error:
46         isl_local_space_free(ls);
47         isl_vec_free(v);
48         return NULL;
49 }
50
51 __isl_give isl_aff *isl_aff_alloc(__isl_take isl_local_space *ls)
52 {
53         isl_ctx *ctx;
54         isl_vec *v;
55         unsigned total;
56
57         if (!ls)
58                 return NULL;
59
60         ctx = isl_local_space_get_ctx(ls);
61         if (!isl_local_space_divs_known(ls))
62                 isl_die(ctx, isl_error_invalid, "local space has unknown divs",
63                         goto error);
64         if (!isl_local_space_is_set(ls))
65                 isl_die(ctx, isl_error_invalid,
66                         "domain of affine expression should be a set",
67                         goto error);
68
69         total = isl_local_space_dim(ls, isl_dim_all);
70         v = isl_vec_alloc(ctx, 1 + 1 + total);
71         return isl_aff_alloc_vec(ls, v);
72 error:
73         isl_local_space_free(ls);
74         return NULL;
75 }
76
77 __isl_give isl_aff *isl_aff_zero_on_domain(__isl_take isl_local_space *ls)
78 {
79         isl_aff *aff;
80
81         aff = isl_aff_alloc(ls);
82         if (!aff)
83                 return NULL;
84
85         isl_int_set_si(aff->v->el[0], 1);
86         isl_seq_clr(aff->v->el + 1, aff->v->size - 1);
87
88         return aff;
89 }
90
91 __isl_give isl_aff *isl_aff_copy(__isl_keep isl_aff *aff)
92 {
93         if (!aff)
94                 return NULL;
95
96         aff->ref++;
97         return aff;
98 }
99
100 __isl_give isl_aff *isl_aff_dup(__isl_keep isl_aff *aff)
101 {
102         if (!aff)
103                 return NULL;
104
105         return isl_aff_alloc_vec(isl_local_space_copy(aff->ls),
106                                  isl_vec_copy(aff->v));
107 }
108
109 __isl_give isl_aff *isl_aff_cow(__isl_take isl_aff *aff)
110 {
111         if (!aff)
112                 return NULL;
113
114         if (aff->ref == 1)
115                 return aff;
116         aff->ref--;
117         return isl_aff_dup(aff);
118 }
119
120 void *isl_aff_free(__isl_take isl_aff *aff)
121 {
122         if (!aff)
123                 return NULL;
124
125         if (--aff->ref > 0)
126                 return NULL;
127
128         isl_local_space_free(aff->ls);
129         isl_vec_free(aff->v);
130
131         free(aff);
132
133         return NULL;
134 }
135
136 isl_ctx *isl_aff_get_ctx(__isl_keep isl_aff *aff)
137 {
138         return aff ? isl_local_space_get_ctx(aff->ls) : NULL;
139 }
140
141 /* Externally, an isl_aff has a map space, but internally, the
142  * ls field corresponds to the domain of that space.
143  */
144 int isl_aff_dim(__isl_keep isl_aff *aff, enum isl_dim_type type)
145 {
146         if (!aff)
147                 return 0;
148         if (type == isl_dim_out)
149                 return 1;
150         if (type == isl_dim_in)
151                 type = isl_dim_set;
152         return isl_local_space_dim(aff->ls, type);
153 }
154
155 __isl_give isl_space *isl_aff_get_domain_space(__isl_keep isl_aff *aff)
156 {
157         return aff ? isl_local_space_get_space(aff->ls) : NULL;
158 }
159
160 __isl_give isl_space *isl_aff_get_space(__isl_keep isl_aff *aff)
161 {
162         isl_space *space;
163         if (!aff)
164                 return NULL;
165         space = isl_local_space_get_space(aff->ls);
166         space = isl_space_from_domain(space);
167         space = isl_space_add_dims(space, isl_dim_out, 1);
168         return space;
169 }
170
171 __isl_give isl_local_space *isl_aff_get_domain_local_space(
172         __isl_keep isl_aff *aff)
173 {
174         return aff ? isl_local_space_copy(aff->ls) : NULL;
175 }
176
177 __isl_give isl_local_space *isl_aff_get_local_space(__isl_keep isl_aff *aff)
178 {
179         isl_local_space *ls;
180         if (!aff)
181                 return NULL;
182         ls = isl_local_space_copy(aff->ls);
183         ls = isl_local_space_from_domain(ls);
184         ls = isl_local_space_add_dims(ls, isl_dim_out, 1);
185         return ls;
186 }
187
188 /* Externally, an isl_aff has a map space, but internally, the
189  * ls field corresponds to the domain of that space.
190  */
191 const char *isl_aff_get_dim_name(__isl_keep isl_aff *aff,
192         enum isl_dim_type type, unsigned pos)
193 {
194         if (!aff)
195                 return NULL;
196         if (type == isl_dim_out)
197                 return NULL;
198         if (type == isl_dim_in)
199                 type = isl_dim_set;
200         return isl_local_space_get_dim_name(aff->ls, type, pos);
201 }
202
203 __isl_give isl_aff *isl_aff_reset_domain_space(__isl_take isl_aff *aff,
204         __isl_take isl_space *dim)
205 {
206         aff = isl_aff_cow(aff);
207         if (!aff || !dim)
208                 goto error;
209
210         aff->ls = isl_local_space_reset_space(aff->ls, dim);
211         if (!aff->ls)
212                 return isl_aff_free(aff);
213
214         return aff;
215 error:
216         isl_aff_free(aff);
217         isl_space_free(dim);
218         return NULL;
219 }
220
221 /* Reset the space of "aff".  This function is called from isl_pw_templ.c
222  * and doesn't know if the space of an element object is represented
223  * directly or through its domain.  It therefore passes along both.
224  */
225 __isl_give isl_aff *isl_aff_reset_space_and_domain(__isl_take isl_aff *aff,
226         __isl_take isl_space *space, __isl_take isl_space *domain)
227 {
228         isl_space_free(space);
229         return isl_aff_reset_domain_space(aff, domain);
230 }
231
232 /* Reorder the coefficients of the affine expression based
233  * on the given reodering.
234  * The reordering r is assumed to have been extended with the local
235  * variables.
236  */
237 static __isl_give isl_vec *vec_reorder(__isl_take isl_vec *vec,
238         __isl_take isl_reordering *r, int n_div)
239 {
240         isl_vec *res;
241         int i;
242
243         if (!vec || !r)
244                 goto error;
245
246         res = isl_vec_alloc(vec->ctx,
247                             2 + isl_space_dim(r->dim, isl_dim_all) + n_div);
248         isl_seq_cpy(res->el, vec->el, 2);
249         isl_seq_clr(res->el + 2, res->size - 2);
250         for (i = 0; i < r->len; ++i)
251                 isl_int_set(res->el[2 + r->pos[i]], vec->el[2 + i]);
252
253         isl_reordering_free(r);
254         isl_vec_free(vec);
255         return res;
256 error:
257         isl_vec_free(vec);
258         isl_reordering_free(r);
259         return NULL;
260 }
261
262 /* Reorder the dimensions of the domain of "aff" according
263  * to the given reordering.
264  */
265 __isl_give isl_aff *isl_aff_realign_domain(__isl_take isl_aff *aff,
266         __isl_take isl_reordering *r)
267 {
268         aff = isl_aff_cow(aff);
269         if (!aff)
270                 goto error;
271
272         r = isl_reordering_extend(r, aff->ls->div->n_row);
273         aff->v = vec_reorder(aff->v, isl_reordering_copy(r),
274                                 aff->ls->div->n_row);
275         aff->ls = isl_local_space_realign(aff->ls, r);
276
277         if (!aff->v || !aff->ls)
278                 return isl_aff_free(aff);
279
280         return aff;
281 error:
282         isl_aff_free(aff);
283         isl_reordering_free(r);
284         return NULL;
285 }
286
287 __isl_give isl_aff *isl_aff_align_params(__isl_take isl_aff *aff,
288         __isl_take isl_space *model)
289 {
290         if (!aff || !model)
291                 goto error;
292
293         if (!isl_space_match(aff->ls->dim, isl_dim_param,
294                              model, isl_dim_param)) {
295                 isl_reordering *exp;
296
297                 model = isl_space_drop_dims(model, isl_dim_in,
298                                         0, isl_space_dim(model, isl_dim_in));
299                 model = isl_space_drop_dims(model, isl_dim_out,
300                                         0, isl_space_dim(model, isl_dim_out));
301                 exp = isl_parameter_alignment_reordering(aff->ls->dim, model);
302                 exp = isl_reordering_extend_space(exp,
303                                         isl_aff_get_domain_space(aff));
304                 aff = isl_aff_realign_domain(aff, exp);
305         }
306
307         isl_space_free(model);
308         return aff;
309 error:
310         isl_space_free(model);
311         isl_aff_free(aff);
312         return NULL;
313 }
314
315 int isl_aff_plain_is_zero(__isl_keep isl_aff *aff)
316 {
317         if (!aff)
318                 return -1;
319
320         return isl_seq_first_non_zero(aff->v->el + 1, aff->v->size - 1) < 0;
321 }
322
323 int isl_aff_plain_is_equal(__isl_keep isl_aff *aff1, __isl_keep isl_aff *aff2)
324 {
325         int equal;
326
327         if (!aff1 || !aff2)
328                 return -1;
329
330         equal = isl_local_space_is_equal(aff1->ls, aff2->ls);
331         if (equal < 0 || !equal)
332                 return equal;
333
334         return isl_vec_is_equal(aff1->v, aff2->v);
335 }
336
337 int isl_aff_get_denominator(__isl_keep isl_aff *aff, isl_int *v)
338 {
339         if (!aff)
340                 return -1;
341         isl_int_set(*v, aff->v->el[0]);
342         return 0;
343 }
344
345 int isl_aff_get_constant(__isl_keep isl_aff *aff, isl_int *v)
346 {
347         if (!aff)
348                 return -1;
349         isl_int_set(*v, aff->v->el[1]);
350         return 0;
351 }
352
353 int isl_aff_get_coefficient(__isl_keep isl_aff *aff,
354         enum isl_dim_type type, int pos, isl_int *v)
355 {
356         if (!aff)
357                 return -1;
358
359         if (type == isl_dim_out)
360                 isl_die(aff->v->ctx, isl_error_invalid,
361                         "output/set dimension does not have a coefficient",
362                         return -1);
363         if (type == isl_dim_in)
364                 type = isl_dim_set;
365
366         if (pos >= isl_local_space_dim(aff->ls, type))
367                 isl_die(aff->v->ctx, isl_error_invalid,
368                         "position out of bounds", return -1);
369
370         pos += isl_local_space_offset(aff->ls, type);
371         isl_int_set(*v, aff->v->el[1 + pos]);
372
373         return 0;
374 }
375
376 __isl_give isl_aff *isl_aff_set_denominator(__isl_take isl_aff *aff, isl_int v)
377 {
378         aff = isl_aff_cow(aff);
379         if (!aff)
380                 return NULL;
381
382         aff->v = isl_vec_cow(aff->v);
383         if (!aff->v)
384                 return isl_aff_free(aff);
385
386         isl_int_set(aff->v->el[0], v);
387
388         return aff;
389 }
390
391 __isl_give isl_aff *isl_aff_set_constant(__isl_take isl_aff *aff, isl_int v)
392 {
393         aff = isl_aff_cow(aff);
394         if (!aff)
395                 return NULL;
396
397         aff->v = isl_vec_cow(aff->v);
398         if (!aff->v)
399                 return isl_aff_free(aff);
400
401         isl_int_set(aff->v->el[1], v);
402
403         return aff;
404 }
405
406 __isl_give isl_aff *isl_aff_add_constant(__isl_take isl_aff *aff, isl_int v)
407 {
408         if (isl_int_is_zero(v))
409                 return aff;
410
411         aff = isl_aff_cow(aff);
412         if (!aff)
413                 return NULL;
414
415         aff->v = isl_vec_cow(aff->v);
416         if (!aff->v)
417                 return isl_aff_free(aff);
418
419         isl_int_addmul(aff->v->el[1], aff->v->el[0], v);
420
421         return aff;
422 }
423
424 __isl_give isl_aff *isl_aff_add_constant_si(__isl_take isl_aff *aff, int v)
425 {
426         isl_int t;
427
428         isl_int_init(t);
429         isl_int_set_si(t, v);
430         aff = isl_aff_add_constant(aff, t);
431         isl_int_clear(t);
432
433         return aff;
434 }
435
436 /* Add "v" to the numerator of the constant term of "aff".
437  */
438 __isl_give isl_aff *isl_aff_add_constant_num(__isl_take isl_aff *aff, isl_int v)
439 {
440         if (isl_int_is_zero(v))
441                 return aff;
442
443         aff = isl_aff_cow(aff);
444         if (!aff)
445                 return NULL;
446
447         aff->v = isl_vec_cow(aff->v);
448         if (!aff->v)
449                 return isl_aff_free(aff);
450
451         isl_int_add(aff->v->el[1], aff->v->el[1], v);
452
453         return aff;
454 }
455
456 /* Add "v" to the numerator of the constant term of "aff".
457  */
458 __isl_give isl_aff *isl_aff_add_constant_num_si(__isl_take isl_aff *aff, int v)
459 {
460         isl_int t;
461
462         if (v == 0)
463                 return aff;
464
465         isl_int_init(t);
466         isl_int_set_si(t, v);
467         aff = isl_aff_add_constant_num(aff, t);
468         isl_int_clear(t);
469
470         return aff;
471 }
472
473 __isl_give isl_aff *isl_aff_set_constant_si(__isl_take isl_aff *aff, int v)
474 {
475         aff = isl_aff_cow(aff);
476         if (!aff)
477                 return NULL;
478
479         aff->v = isl_vec_cow(aff->v);
480         if (!aff->v)
481                 return isl_aff_free(aff);
482
483         isl_int_set_si(aff->v->el[1], v);
484
485         return aff;
486 }
487
488 __isl_give isl_aff *isl_aff_set_coefficient(__isl_take isl_aff *aff,
489         enum isl_dim_type type, int pos, isl_int v)
490 {
491         if (!aff)
492                 return NULL;
493
494         if (type == isl_dim_out)
495                 isl_die(aff->v->ctx, isl_error_invalid,
496                         "output/set dimension does not have a coefficient",
497                         return isl_aff_free(aff));
498         if (type == isl_dim_in)
499                 type = isl_dim_set;
500
501         if (pos >= isl_local_space_dim(aff->ls, type))
502                 isl_die(aff->v->ctx, isl_error_invalid,
503                         "position out of bounds", return isl_aff_free(aff));
504
505         aff = isl_aff_cow(aff);
506         if (!aff)
507                 return NULL;
508
509         aff->v = isl_vec_cow(aff->v);
510         if (!aff->v)
511                 return isl_aff_free(aff);
512
513         pos += isl_local_space_offset(aff->ls, type);
514         isl_int_set(aff->v->el[1 + pos], v);
515
516         return aff;
517 }
518
519 __isl_give isl_aff *isl_aff_set_coefficient_si(__isl_take isl_aff *aff,
520         enum isl_dim_type type, int pos, int v)
521 {
522         if (!aff)
523                 return NULL;
524
525         if (type == isl_dim_out)
526                 isl_die(aff->v->ctx, isl_error_invalid,
527                         "output/set dimension does not have a coefficient",
528                         return isl_aff_free(aff));
529         if (type == isl_dim_in)
530                 type = isl_dim_set;
531
532         if (pos >= isl_local_space_dim(aff->ls, type))
533                 isl_die(aff->v->ctx, isl_error_invalid,
534                         "position out of bounds", return isl_aff_free(aff));
535
536         aff = isl_aff_cow(aff);
537         if (!aff)
538                 return NULL;
539
540         aff->v = isl_vec_cow(aff->v);
541         if (!aff->v)
542                 return isl_aff_free(aff);
543
544         pos += isl_local_space_offset(aff->ls, type);
545         isl_int_set_si(aff->v->el[1 + pos], v);
546
547         return aff;
548 }
549
550 __isl_give isl_aff *isl_aff_add_coefficient(__isl_take isl_aff *aff,
551         enum isl_dim_type type, int pos, isl_int v)
552 {
553         if (!aff)
554                 return NULL;
555
556         if (type == isl_dim_out)
557                 isl_die(aff->v->ctx, isl_error_invalid,
558                         "output/set dimension does not have a coefficient",
559                         return isl_aff_free(aff));
560         if (type == isl_dim_in)
561                 type = isl_dim_set;
562
563         if (pos >= isl_local_space_dim(aff->ls, type))
564                 isl_die(aff->v->ctx, isl_error_invalid,
565                         "position out of bounds", return isl_aff_free(aff));
566
567         aff = isl_aff_cow(aff);
568         if (!aff)
569                 return NULL;
570
571         aff->v = isl_vec_cow(aff->v);
572         if (!aff->v)
573                 return isl_aff_free(aff);
574
575         pos += isl_local_space_offset(aff->ls, type);
576         isl_int_addmul(aff->v->el[1 + pos], aff->v->el[0], v);
577
578         return aff;
579 }
580
581 __isl_give isl_aff *isl_aff_add_coefficient_si(__isl_take isl_aff *aff,
582         enum isl_dim_type type, int pos, int v)
583 {
584         isl_int t;
585
586         isl_int_init(t);
587         isl_int_set_si(t, v);
588         aff = isl_aff_add_coefficient(aff, type, pos, t);
589         isl_int_clear(t);
590
591         return aff;
592 }
593
594 __isl_give isl_aff *isl_aff_get_div(__isl_keep isl_aff *aff, int pos)
595 {
596         if (!aff)
597                 return NULL;
598
599         return isl_local_space_get_div(aff->ls, pos);
600 }
601
602 __isl_give isl_aff *isl_aff_neg(__isl_take isl_aff *aff)
603 {
604         aff = isl_aff_cow(aff);
605         if (!aff)
606                 return NULL;
607         aff->v = isl_vec_cow(aff->v);
608         if (!aff->v)
609                 return isl_aff_free(aff);
610
611         isl_seq_neg(aff->v->el + 1, aff->v->el + 1, aff->v->size - 1);
612
613         return aff;
614 }
615
616 /* Remove divs from the local space that do not appear in the affine
617  * expression.
618  * We currently only remove divs at the end.
619  * Some intermediate divs may also not appear directly in the affine
620  * expression, but we would also need to check that no other divs are
621  * defined in terms of them.
622  */
623 __isl_give isl_aff *isl_aff_remove_unused_divs( __isl_take isl_aff *aff)
624 {
625         int pos;
626         int off;
627         int n;
628
629         if (!aff)
630                 return NULL;
631
632         n = isl_local_space_dim(aff->ls, isl_dim_div);
633         off = isl_local_space_offset(aff->ls, isl_dim_div);
634
635         pos = isl_seq_last_non_zero(aff->v->el + 1 + off, n) + 1;
636         if (pos == n)
637                 return aff;
638
639         aff = isl_aff_cow(aff);
640         if (!aff)
641                 return NULL;
642
643         aff->ls = isl_local_space_drop_dims(aff->ls, isl_dim_div, pos, n - pos);
644         aff->v = isl_vec_drop_els(aff->v, 1 + off + pos, n - pos);
645         if (!aff->ls || !aff->v)
646                 return isl_aff_free(aff);
647
648         return aff;
649 }
650
651 __isl_give isl_aff *isl_aff_normalize(__isl_take isl_aff *aff)
652 {
653         if (!aff)
654                 return NULL;
655         aff->v = isl_vec_normalize(aff->v);
656         if (!aff->v)
657                 return isl_aff_free(aff);
658         aff = isl_aff_remove_unused_divs(aff);
659         return aff;
660 }
661
662 /* Given f, return floor(f).
663  * If f is an integer expression, then just return f.
664  * If f is a constant, then return the constant floor(f).
665  * Otherwise, if f = g/m, write g = q m + r,
666  * create a new div d = [r/m] and return the expression q + d.
667  * The coefficients in r are taken to lie between -m/2 and m/2.
668  */
669 __isl_give isl_aff *isl_aff_floor(__isl_take isl_aff *aff)
670 {
671         int i;
672         int size;
673         isl_ctx *ctx;
674         isl_vec *div;
675
676         if (!aff)
677                 return NULL;
678
679         if (isl_int_is_one(aff->v->el[0]))
680                 return aff;
681
682         aff = isl_aff_cow(aff);
683         if (!aff)
684                 return NULL;
685
686         aff->v = isl_vec_cow(aff->v);
687         if (!aff->v)
688                 return isl_aff_free(aff);
689
690         if (isl_aff_is_cst(aff)) {
691                 isl_int_fdiv_q(aff->v->el[1], aff->v->el[1], aff->v->el[0]);
692                 isl_int_set_si(aff->v->el[0], 1);
693                 return aff;
694         }
695
696         div = isl_vec_copy(aff->v);
697         div = isl_vec_cow(div);
698         if (!div)
699                 return isl_aff_free(aff);
700
701         ctx = isl_aff_get_ctx(aff);
702         isl_int_fdiv_q(aff->v->el[0], aff->v->el[0], ctx->two);
703         for (i = 1; i < aff->v->size; ++i) {
704                 isl_int_fdiv_r(div->el[i], div->el[i], div->el[0]);
705                 isl_int_fdiv_q(aff->v->el[i], aff->v->el[i], div->el[0]);
706                 if (isl_int_gt(div->el[i], aff->v->el[0])) {
707                         isl_int_sub(div->el[i], div->el[i], div->el[0]);
708                         isl_int_add_ui(aff->v->el[i], aff->v->el[i], 1);
709                 }
710         }
711
712         aff->ls = isl_local_space_add_div(aff->ls, div);
713         if (!aff->ls)
714                 return isl_aff_free(aff);
715
716         size = aff->v->size;
717         aff->v = isl_vec_extend(aff->v, size + 1);
718         if (!aff->v)
719                 return isl_aff_free(aff);
720         isl_int_set_si(aff->v->el[0], 1);
721         isl_int_set_si(aff->v->el[size], 1);
722
723         return aff;
724 }
725
726 /* Compute
727  *
728  *      aff mod m = aff - m * floor(aff/m)
729  */
730 __isl_give isl_aff *isl_aff_mod(__isl_take isl_aff *aff, isl_int m)
731 {
732         isl_aff *res;
733
734         res = isl_aff_copy(aff);
735         aff = isl_aff_scale_down(aff, m);
736         aff = isl_aff_floor(aff);
737         aff = isl_aff_scale(aff, m);
738         res = isl_aff_sub(res, aff);
739
740         return res;
741 }
742
743 /* Compute
744  *
745  *      pwaff mod m = pwaff - m * floor(pwaff/m)
746  */
747 __isl_give isl_pw_aff *isl_pw_aff_mod(__isl_take isl_pw_aff *pwaff, isl_int m)
748 {
749         isl_pw_aff *res;
750
751         res = isl_pw_aff_copy(pwaff);
752         pwaff = isl_pw_aff_scale_down(pwaff, m);
753         pwaff = isl_pw_aff_floor(pwaff);
754         pwaff = isl_pw_aff_scale(pwaff, m);
755         res = isl_pw_aff_sub(res, pwaff);
756
757         return res;
758 }
759
760 /* Given f, return ceil(f).
761  * If f is an integer expression, then just return f.
762  * Otherwise, create a new div d = [-f] and return the expression -d.
763  */
764 __isl_give isl_aff *isl_aff_ceil(__isl_take isl_aff *aff)
765 {
766         if (!aff)
767                 return NULL;
768
769         if (isl_int_is_one(aff->v->el[0]))
770                 return aff;
771
772         aff = isl_aff_neg(aff);
773         aff = isl_aff_floor(aff);
774         aff = isl_aff_neg(aff);
775
776         return aff;
777 }
778
779 /* Apply the expansion computed by isl_merge_divs.
780  * The expansion itself is given by "exp" while the resulting
781  * list of divs is given by "div".
782  */
783 __isl_give isl_aff *isl_aff_expand_divs( __isl_take isl_aff *aff,
784         __isl_take isl_mat *div, int *exp)
785 {
786         int i, j;
787         int old_n_div;
788         int new_n_div;
789         int offset;
790
791         aff = isl_aff_cow(aff);
792         if (!aff || !div)
793                 goto error;
794
795         old_n_div = isl_local_space_dim(aff->ls, isl_dim_div);
796         new_n_div = isl_mat_rows(div);
797         if (new_n_div < old_n_div)
798                 isl_die(isl_mat_get_ctx(div), isl_error_invalid,
799                         "not an expansion", goto error);
800
801         aff->v = isl_vec_extend(aff->v, aff->v->size + new_n_div - old_n_div);
802         if (!aff->v)
803                 goto error;
804
805         offset = 1 + isl_local_space_offset(aff->ls, isl_dim_div);
806         j = old_n_div - 1;
807         for (i = new_n_div - 1; i >= 0; --i) {
808                 if (j >= 0 && exp[j] == i) {
809                         if (i != j)
810                                 isl_int_swap(aff->v->el[offset + i],
811                                              aff->v->el[offset + j]);
812                         j--;
813                 } else
814                         isl_int_set_si(aff->v->el[offset + i], 0);
815         }
816
817         aff->ls = isl_local_space_replace_divs(aff->ls, isl_mat_copy(div));
818         if (!aff->ls)
819                 goto error;
820         isl_mat_free(div);
821         return aff;
822 error:
823         isl_aff_free(aff);
824         isl_mat_free(div);
825         return NULL;
826 }
827
828 /* Add two affine expressions that live in the same local space.
829  */
830 static __isl_give isl_aff *add_expanded(__isl_take isl_aff *aff1,
831         __isl_take isl_aff *aff2)
832 {
833         isl_int gcd, f;
834
835         aff1 = isl_aff_cow(aff1);
836         if (!aff1 || !aff2)
837                 goto error;
838
839         aff1->v = isl_vec_cow(aff1->v);
840         if (!aff1->v)
841                 goto error;
842
843         isl_int_init(gcd);
844         isl_int_init(f);
845         isl_int_gcd(gcd, aff1->v->el[0], aff2->v->el[0]);
846         isl_int_divexact(f, aff2->v->el[0], gcd);
847         isl_seq_scale(aff1->v->el + 1, aff1->v->el + 1, f, aff1->v->size - 1);
848         isl_int_divexact(f, aff1->v->el[0], gcd);
849         isl_seq_addmul(aff1->v->el + 1, f, aff2->v->el + 1, aff1->v->size - 1);
850         isl_int_divexact(f, aff2->v->el[0], gcd);
851         isl_int_mul(aff1->v->el[0], aff1->v->el[0], f);
852         isl_int_clear(f);
853         isl_int_clear(gcd);
854
855         isl_aff_free(aff2);
856         return aff1;
857 error:
858         isl_aff_free(aff1);
859         isl_aff_free(aff2);
860         return NULL;
861 }
862
863 __isl_give isl_aff *isl_aff_add(__isl_take isl_aff *aff1,
864         __isl_take isl_aff *aff2)
865 {
866         isl_ctx *ctx;
867         int *exp1 = NULL;
868         int *exp2 = NULL;
869         isl_mat *div;
870
871         if (!aff1 || !aff2)
872                 goto error;
873
874         ctx = isl_aff_get_ctx(aff1);
875         if (!isl_space_is_equal(aff1->ls->dim, aff2->ls->dim))
876                 isl_die(ctx, isl_error_invalid,
877                         "spaces don't match", goto error);
878
879         if (aff1->ls->div->n_row == 0 && aff2->ls->div->n_row == 0)
880                 return add_expanded(aff1, aff2);
881
882         exp1 = isl_alloc_array(ctx, int, aff1->ls->div->n_row);
883         exp2 = isl_alloc_array(ctx, int, aff2->ls->div->n_row);
884         if (!exp1 || !exp2)
885                 goto error;
886
887         div = isl_merge_divs(aff1->ls->div, aff2->ls->div, exp1, exp2);
888         aff1 = isl_aff_expand_divs(aff1, isl_mat_copy(div), exp1);
889         aff2 = isl_aff_expand_divs(aff2, div, exp2);
890         free(exp1);
891         free(exp2);
892
893         return add_expanded(aff1, aff2);
894 error:
895         free(exp1);
896         free(exp2);
897         isl_aff_free(aff1);
898         isl_aff_free(aff2);
899         return NULL;
900 }
901
902 __isl_give isl_aff *isl_aff_sub(__isl_take isl_aff *aff1,
903         __isl_take isl_aff *aff2)
904 {
905         return isl_aff_add(aff1, isl_aff_neg(aff2));
906 }
907
908 __isl_give isl_aff *isl_aff_scale(__isl_take isl_aff *aff, isl_int f)
909 {
910         isl_int gcd;
911
912         if (isl_int_is_one(f))
913                 return aff;
914
915         aff = isl_aff_cow(aff);
916         if (!aff)
917                 return NULL;
918         aff->v = isl_vec_cow(aff->v);
919         if (!aff->v)
920                 return isl_aff_free(aff);
921
922         isl_int_init(gcd);
923         isl_int_gcd(gcd, aff->v->el[0], f);
924         isl_int_divexact(aff->v->el[0], aff->v->el[0], gcd);
925         isl_int_divexact(gcd, f, gcd);
926         isl_seq_scale(aff->v->el + 1, aff->v->el + 1, gcd, aff->v->size - 1);
927         isl_int_clear(gcd);
928
929         return aff;
930 }
931
932 __isl_give isl_aff *isl_aff_scale_down(__isl_take isl_aff *aff, isl_int f)
933 {
934         isl_int gcd;
935
936         if (isl_int_is_one(f))
937                 return aff;
938
939         aff = isl_aff_cow(aff);
940         if (!aff)
941                 return NULL;
942         aff->v = isl_vec_cow(aff->v);
943         if (!aff->v)
944                 return isl_aff_free(aff);
945
946         isl_int_init(gcd);
947         isl_seq_gcd(aff->v->el + 1, aff->v->size - 1, &gcd);
948         isl_int_gcd(gcd, gcd, f);
949         isl_seq_scale_down(aff->v->el + 1, aff->v->el + 1, gcd, aff->v->size - 1);
950         isl_int_divexact(gcd, f, gcd);
951         isl_int_mul(aff->v->el[0], aff->v->el[0], gcd);
952         isl_int_clear(gcd);
953
954         return aff;
955 }
956
957 __isl_give isl_aff *isl_aff_scale_down_ui(__isl_take isl_aff *aff, unsigned f)
958 {
959         isl_int v;
960
961         if (f == 1)
962                 return aff;
963
964         isl_int_init(v);
965         isl_int_set_ui(v, f);
966         aff = isl_aff_scale_down(aff, v);
967         isl_int_clear(v);
968
969         return aff;
970 }
971
972 __isl_give isl_aff *isl_aff_set_dim_name(__isl_take isl_aff *aff,
973         enum isl_dim_type type, unsigned pos, const char *s)
974 {
975         aff = isl_aff_cow(aff);
976         if (!aff)
977                 return NULL;
978         if (type == isl_dim_out)
979                 isl_die(aff->v->ctx, isl_error_invalid,
980                         "cannot set name of output/set dimension",
981                         return isl_aff_free(aff));
982         if (type == isl_dim_in)
983                 type = isl_dim_set;
984         aff->ls = isl_local_space_set_dim_name(aff->ls, type, pos, s);
985         if (!aff->ls)
986                 return isl_aff_free(aff);
987
988         return aff;
989 }
990
991 __isl_give isl_aff *isl_aff_set_dim_id(__isl_take isl_aff *aff,
992         enum isl_dim_type type, unsigned pos, __isl_take isl_id *id)
993 {
994         aff = isl_aff_cow(aff);
995         if (!aff)
996                 return isl_id_free(id);
997         if (type == isl_dim_out)
998                 isl_die(aff->v->ctx, isl_error_invalid,
999                         "cannot set name of output/set dimension",
1000                         goto error);
1001         if (type == isl_dim_in)
1002                 type = isl_dim_set;
1003         aff->ls = isl_local_space_set_dim_id(aff->ls, type, pos, id);
1004         if (!aff->ls)
1005                 return isl_aff_free(aff);
1006
1007         return aff;
1008 error:
1009         isl_id_free(id);
1010         isl_aff_free(aff);
1011         return NULL;
1012 }
1013
1014 /* Exploit the equalities in "eq" to simplify the affine expression
1015  * and the expressions of the integer divisions in the local space.
1016  * The integer divisions in this local space are assumed to appear
1017  * as regular dimensions in "eq".
1018  */
1019 static __isl_give isl_aff *isl_aff_substitute_equalities_lifted(
1020         __isl_take isl_aff *aff, __isl_take isl_basic_set *eq)
1021 {
1022         int i, j;
1023         unsigned total;
1024         unsigned n_div;
1025
1026         if (!eq)
1027                 goto error;
1028         if (eq->n_eq == 0) {
1029                 isl_basic_set_free(eq);
1030                 return aff;
1031         }
1032
1033         aff = isl_aff_cow(aff);
1034         if (!aff)
1035                 goto error;
1036
1037         aff->ls = isl_local_space_substitute_equalities(aff->ls,
1038                                                         isl_basic_set_copy(eq));
1039         if (!aff->ls)
1040                 goto error;
1041
1042         total = 1 + isl_space_dim(eq->dim, isl_dim_all);
1043         n_div = eq->n_div;
1044         for (i = 0; i < eq->n_eq; ++i) {
1045                 j = isl_seq_last_non_zero(eq->eq[i], total + n_div);
1046                 if (j < 0 || j == 0 || j >= total)
1047                         continue;
1048
1049                 isl_seq_elim(aff->v->el + 1, eq->eq[i], j, total,
1050                                 &aff->v->el[0]);
1051         }
1052
1053         isl_basic_set_free(eq);
1054         aff = isl_aff_normalize(aff);
1055         return aff;
1056 error:
1057         isl_basic_set_free(eq);
1058         isl_aff_free(aff);
1059         return NULL;
1060 }
1061
1062 /* Exploit the equalities in "eq" to simplify the affine expression
1063  * and the expressions of the integer divisions in the local space.
1064  */
1065 static __isl_give isl_aff *isl_aff_substitute_equalities(
1066         __isl_take isl_aff *aff, __isl_take isl_basic_set *eq)
1067 {
1068         int n_div;
1069
1070         if (!aff || !eq)
1071                 goto error;
1072         n_div = isl_local_space_dim(aff->ls, isl_dim_div);
1073         if (n_div > 0)
1074                 eq = isl_basic_set_add(eq, isl_dim_set, n_div);
1075         return isl_aff_substitute_equalities_lifted(aff, eq);
1076 error:
1077         isl_basic_set_free(eq);
1078         isl_aff_free(aff);
1079         return NULL;
1080 }
1081
1082 /* Look for equalities among the variables shared by context and aff
1083  * and the integer divisions of aff, if any.
1084  * The equalities are then used to eliminate coefficients and/or integer
1085  * divisions from aff.
1086  */
1087 __isl_give isl_aff *isl_aff_gist(__isl_take isl_aff *aff,
1088         __isl_take isl_set *context)
1089 {
1090         isl_basic_set *hull;
1091         int n_div;
1092
1093         if (!aff)
1094                 goto error;
1095         n_div = isl_local_space_dim(aff->ls, isl_dim_div);
1096         if (n_div > 0) {
1097                 isl_basic_set *bset;
1098                 isl_local_space *ls;
1099                 context = isl_set_add_dims(context, isl_dim_set, n_div);
1100                 ls = isl_aff_get_domain_local_space(aff);
1101                 bset = isl_basic_set_from_local_space(ls);
1102                 bset = isl_basic_set_lift(bset);
1103                 bset = isl_basic_set_flatten(bset);
1104                 context = isl_set_intersect(context,
1105                                             isl_set_from_basic_set(bset));
1106         }
1107
1108         hull = isl_set_affine_hull(context);
1109         return isl_aff_substitute_equalities_lifted(aff, hull);
1110 error:
1111         isl_aff_free(aff);
1112         isl_set_free(context);
1113         return NULL;
1114 }
1115
1116 __isl_give isl_aff *isl_aff_gist_params(__isl_take isl_aff *aff,
1117         __isl_take isl_set *context)
1118 {
1119         isl_set *dom_context = isl_set_universe(isl_aff_get_domain_space(aff));
1120         dom_context = isl_set_intersect_params(dom_context, context);
1121         return isl_aff_gist(aff, dom_context);
1122 }
1123
1124 /* Return a basic set containing those elements in the space
1125  * of aff where it is non-negative.
1126  */
1127 __isl_give isl_basic_set *isl_aff_nonneg_basic_set(__isl_take isl_aff *aff)
1128 {
1129         isl_constraint *ineq;
1130         isl_basic_set *bset;
1131
1132         ineq = isl_inequality_from_aff(aff);
1133
1134         bset = isl_basic_set_from_constraint(ineq);
1135         bset = isl_basic_set_simplify(bset);
1136         return bset;
1137 }
1138
1139 /* Return a basic set containing those elements in the space
1140  * of aff where it is zero.
1141  */
1142 __isl_give isl_basic_set *isl_aff_zero_basic_set(__isl_take isl_aff *aff)
1143 {
1144         isl_constraint *ineq;
1145         isl_basic_set *bset;
1146
1147         ineq = isl_equality_from_aff(aff);
1148
1149         bset = isl_basic_set_from_constraint(ineq);
1150         bset = isl_basic_set_simplify(bset);
1151         return bset;
1152 }
1153
1154 /* Return a basic set containing those elements in the shared space
1155  * of aff1 and aff2 where aff1 is greater than or equal to aff2.
1156  */
1157 __isl_give isl_basic_set *isl_aff_ge_basic_set(__isl_take isl_aff *aff1,
1158         __isl_take isl_aff *aff2)
1159 {
1160         aff1 = isl_aff_sub(aff1, aff2);
1161
1162         return isl_aff_nonneg_basic_set(aff1);
1163 }
1164
1165 /* Return a basic set containing those elements in the shared space
1166  * of aff1 and aff2 where aff1 is smaller than or equal to aff2.
1167  */
1168 __isl_give isl_basic_set *isl_aff_le_basic_set(__isl_take isl_aff *aff1,
1169         __isl_take isl_aff *aff2)
1170 {
1171         return isl_aff_ge_basic_set(aff2, aff1);
1172 }
1173
1174 __isl_give isl_aff *isl_aff_add_on_domain(__isl_keep isl_set *dom,
1175         __isl_take isl_aff *aff1, __isl_take isl_aff *aff2)
1176 {
1177         aff1 = isl_aff_add(aff1, aff2);
1178         aff1 = isl_aff_gist(aff1, isl_set_copy(dom));
1179         return aff1;
1180 }
1181
1182 int isl_aff_is_empty(__isl_keep isl_aff *aff)
1183 {
1184         if (!aff)
1185                 return -1;
1186
1187         return 0;
1188 }
1189
1190 /* Check whether the given affine expression has non-zero coefficient
1191  * for any dimension in the given range or if any of these dimensions
1192  * appear with non-zero coefficients in any of the integer divisions
1193  * involved in the affine expression.
1194  */
1195 int isl_aff_involves_dims(__isl_keep isl_aff *aff,
1196         enum isl_dim_type type, unsigned first, unsigned n)
1197 {
1198         int i;
1199         isl_ctx *ctx;
1200         int *active = NULL;
1201         int involves = 0;
1202
1203         if (!aff)
1204                 return -1;
1205         if (n == 0)
1206                 return 0;
1207
1208         ctx = isl_aff_get_ctx(aff);
1209         if (first + n > isl_aff_dim(aff, type))
1210                 isl_die(ctx, isl_error_invalid,
1211                         "range out of bounds", return -1);
1212
1213         active = isl_local_space_get_active(aff->ls, aff->v->el + 2);
1214         if (!active)
1215                 goto error;
1216
1217         first += isl_local_space_offset(aff->ls, type) - 1;
1218         for (i = 0; i < n; ++i)
1219                 if (active[first + i]) {
1220                         involves = 1;
1221                         break;
1222                 }
1223
1224         free(active);
1225
1226         return involves;
1227 error:
1228         free(active);
1229         return -1;
1230 }
1231
1232 __isl_give isl_aff *isl_aff_drop_dims(__isl_take isl_aff *aff,
1233         enum isl_dim_type type, unsigned first, unsigned n)
1234 {
1235         isl_ctx *ctx;
1236
1237         if (!aff)
1238                 return NULL;
1239         if (type == isl_dim_out)
1240                 isl_die(aff->v->ctx, isl_error_invalid,
1241                         "cannot drop output/set dimension",
1242                         return isl_aff_free(aff));
1243         if (type == isl_dim_in)
1244                 type = isl_dim_set;
1245         if (n == 0 && !isl_local_space_is_named_or_nested(aff->ls, type))
1246                 return aff;
1247
1248         ctx = isl_aff_get_ctx(aff);
1249         if (first + n > isl_local_space_dim(aff->ls, type))
1250                 isl_die(ctx, isl_error_invalid, "range out of bounds",
1251                         return isl_aff_free(aff));
1252
1253         aff = isl_aff_cow(aff);
1254         if (!aff)
1255                 return NULL;
1256
1257         aff->ls = isl_local_space_drop_dims(aff->ls, type, first, n);
1258         if (!aff->ls)
1259                 return isl_aff_free(aff);
1260
1261         first += 1 + isl_local_space_offset(aff->ls, type);
1262         aff->v = isl_vec_drop_els(aff->v, first, n);
1263         if (!aff->v)
1264                 return isl_aff_free(aff);
1265
1266         return aff;
1267 }
1268
1269 /* Project the domain of the affine expression onto its parameter space.
1270  * The affine expression may not involve any of the domain dimensions.
1271  */
1272 __isl_give isl_aff *isl_aff_project_domain_on_params(__isl_take isl_aff *aff)
1273 {
1274         isl_space *space;
1275         unsigned n;
1276         int involves;
1277
1278         n = isl_aff_dim(aff, isl_dim_in);
1279         involves = isl_aff_involves_dims(aff, isl_dim_in, 0, n);
1280         if (involves < 0)
1281                 return isl_aff_free(aff);
1282         if (involves)
1283                 isl_die(isl_aff_get_ctx(aff), isl_error_invalid,
1284                     "affine expression involves some of the domain dimensions",
1285                     return isl_aff_free(aff));
1286         aff = isl_aff_drop_dims(aff, isl_dim_in, 0, n);
1287         space = isl_aff_get_domain_space(aff);
1288         space = isl_space_params(space);
1289         aff = isl_aff_reset_domain_space(aff, space);
1290         return aff;
1291 }
1292
1293 __isl_give isl_aff *isl_aff_insert_dims(__isl_take isl_aff *aff,
1294         enum isl_dim_type type, unsigned first, unsigned n)
1295 {
1296         isl_ctx *ctx;
1297
1298         if (!aff)
1299                 return NULL;
1300         if (type == isl_dim_out)
1301                 isl_die(aff->v->ctx, isl_error_invalid,
1302                         "cannot insert output/set dimensions",
1303                         return isl_aff_free(aff));
1304         if (type == isl_dim_in)
1305                 type = isl_dim_set;
1306         if (n == 0 && !isl_local_space_is_named_or_nested(aff->ls, type))
1307                 return aff;
1308
1309         ctx = isl_aff_get_ctx(aff);
1310         if (first > isl_local_space_dim(aff->ls, type))
1311                 isl_die(ctx, isl_error_invalid, "position out of bounds",
1312                         return isl_aff_free(aff));
1313
1314         aff = isl_aff_cow(aff);
1315         if (!aff)
1316                 return NULL;
1317
1318         aff->ls = isl_local_space_insert_dims(aff->ls, type, first, n);
1319         if (!aff->ls)
1320                 return isl_aff_free(aff);
1321
1322         first += 1 + isl_local_space_offset(aff->ls, type);
1323         aff->v = isl_vec_insert_zero_els(aff->v, first, n);
1324         if (!aff->v)
1325                 return isl_aff_free(aff);
1326
1327         return aff;
1328 }
1329
1330 __isl_give isl_aff *isl_aff_add_dims(__isl_take isl_aff *aff,
1331         enum isl_dim_type type, unsigned n)
1332 {
1333         unsigned pos;
1334
1335         pos = isl_aff_dim(aff, type);
1336
1337         return isl_aff_insert_dims(aff, type, pos, n);
1338 }
1339
1340 __isl_give isl_pw_aff *isl_pw_aff_add_dims(__isl_take isl_pw_aff *pwaff,
1341         enum isl_dim_type type, unsigned n)
1342 {
1343         unsigned pos;
1344
1345         pos = isl_pw_aff_dim(pwaff, type);
1346
1347         return isl_pw_aff_insert_dims(pwaff, type, pos, n);
1348 }
1349
1350 __isl_give isl_pw_aff *isl_pw_aff_from_aff(__isl_take isl_aff *aff)
1351 {
1352         isl_set *dom = isl_set_universe(isl_aff_get_domain_space(aff));
1353         return isl_pw_aff_alloc(dom, aff);
1354 }
1355
1356 #undef PW
1357 #define PW isl_pw_aff
1358 #undef EL
1359 #define EL isl_aff
1360 #undef EL_IS_ZERO
1361 #define EL_IS_ZERO is_empty
1362 #undef ZERO
1363 #define ZERO empty
1364 #undef IS_ZERO
1365 #define IS_ZERO is_empty
1366 #undef FIELD
1367 #define FIELD aff
1368 #undef DEFAULT_IS_ZERO
1369 #define DEFAULT_IS_ZERO 0
1370
1371 #define NO_EVAL
1372 #define NO_OPT
1373 #define NO_MOVE_DIMS
1374 #define NO_LIFT
1375 #define NO_MORPH
1376
1377 #include <isl_pw_templ.c>
1378
1379 static __isl_give isl_set *align_params_pw_pw_set_and(
1380         __isl_take isl_pw_aff *pwaff1, __isl_take isl_pw_aff *pwaff2,
1381         __isl_give isl_set *(*fn)(__isl_take isl_pw_aff *pwaff1,
1382                                     __isl_take isl_pw_aff *pwaff2))
1383 {
1384         if (!pwaff1 || !pwaff2)
1385                 goto error;
1386         if (isl_space_match(pwaff1->dim, isl_dim_param,
1387                           pwaff2->dim, isl_dim_param))
1388                 return fn(pwaff1, pwaff2);
1389         if (!isl_space_has_named_params(pwaff1->dim) ||
1390             !isl_space_has_named_params(pwaff2->dim))
1391                 isl_die(isl_pw_aff_get_ctx(pwaff1), isl_error_invalid,
1392                         "unaligned unnamed parameters", goto error);
1393         pwaff1 = isl_pw_aff_align_params(pwaff1, isl_pw_aff_get_space(pwaff2));
1394         pwaff2 = isl_pw_aff_align_params(pwaff2, isl_pw_aff_get_space(pwaff1));
1395         return fn(pwaff1, pwaff2);
1396 error:
1397         isl_pw_aff_free(pwaff1);
1398         isl_pw_aff_free(pwaff2);
1399         return NULL;
1400 }
1401
1402 /* Compute a piecewise quasi-affine expression with a domain that
1403  * is the union of those of pwaff1 and pwaff2 and such that on each
1404  * cell, the quasi-affine expression is the better (according to cmp)
1405  * of those of pwaff1 and pwaff2.  If only one of pwaff1 or pwaff2
1406  * is defined on a given cell, then the associated expression
1407  * is the defined one.
1408  */
1409 static __isl_give isl_pw_aff *pw_aff_union_opt(__isl_take isl_pw_aff *pwaff1,
1410         __isl_take isl_pw_aff *pwaff2,
1411         __isl_give isl_basic_set *(*cmp)(__isl_take isl_aff *aff1,
1412                                         __isl_take isl_aff *aff2))
1413 {
1414         int i, j, n;
1415         isl_pw_aff *res;
1416         isl_ctx *ctx;
1417         isl_set *set;
1418
1419         if (!pwaff1 || !pwaff2)
1420                 goto error;
1421
1422         ctx = isl_space_get_ctx(pwaff1->dim);
1423         if (!isl_space_is_equal(pwaff1->dim, pwaff2->dim))
1424                 isl_die(ctx, isl_error_invalid,
1425                         "arguments should live in same space", goto error);
1426
1427         if (isl_pw_aff_is_empty(pwaff1)) {
1428                 isl_pw_aff_free(pwaff1);
1429                 return pwaff2;
1430         }
1431
1432         if (isl_pw_aff_is_empty(pwaff2)) {
1433                 isl_pw_aff_free(pwaff2);
1434                 return pwaff1;
1435         }
1436
1437         n = 2 * (pwaff1->n + 1) * (pwaff2->n + 1);
1438         res = isl_pw_aff_alloc_size(isl_space_copy(pwaff1->dim), n);
1439
1440         for (i = 0; i < pwaff1->n; ++i) {
1441                 set = isl_set_copy(pwaff1->p[i].set);
1442                 for (j = 0; j < pwaff2->n; ++j) {
1443                         struct isl_set *common;
1444                         isl_set *better;
1445
1446                         common = isl_set_intersect(
1447                                         isl_set_copy(pwaff1->p[i].set),
1448                                         isl_set_copy(pwaff2->p[j].set));
1449                         better = isl_set_from_basic_set(cmp(
1450                                         isl_aff_copy(pwaff2->p[j].aff),
1451                                         isl_aff_copy(pwaff1->p[i].aff)));
1452                         better = isl_set_intersect(common, better);
1453                         if (isl_set_plain_is_empty(better)) {
1454                                 isl_set_free(better);
1455                                 continue;
1456                         }
1457                         set = isl_set_subtract(set, isl_set_copy(better));
1458
1459                         res = isl_pw_aff_add_piece(res, better,
1460                                                 isl_aff_copy(pwaff2->p[j].aff));
1461                 }
1462                 res = isl_pw_aff_add_piece(res, set,
1463                                                 isl_aff_copy(pwaff1->p[i].aff));
1464         }
1465
1466         for (j = 0; j < pwaff2->n; ++j) {
1467                 set = isl_set_copy(pwaff2->p[j].set);
1468                 for (i = 0; i < pwaff1->n; ++i)
1469                         set = isl_set_subtract(set,
1470                                         isl_set_copy(pwaff1->p[i].set));
1471                 res = isl_pw_aff_add_piece(res, set,
1472                                                 isl_aff_copy(pwaff2->p[j].aff));
1473         }
1474
1475         isl_pw_aff_free(pwaff1);
1476         isl_pw_aff_free(pwaff2);
1477
1478         return res;
1479 error:
1480         isl_pw_aff_free(pwaff1);
1481         isl_pw_aff_free(pwaff2);
1482         return NULL;
1483 }
1484
1485 /* Compute a piecewise quasi-affine expression with a domain that
1486  * is the union of those of pwaff1 and pwaff2 and such that on each
1487  * cell, the quasi-affine expression is the maximum of those of pwaff1
1488  * and pwaff2.  If only one of pwaff1 or pwaff2 is defined on a given
1489  * cell, then the associated expression is the defined one.
1490  */
1491 static __isl_give isl_pw_aff *pw_aff_union_max(__isl_take isl_pw_aff *pwaff1,
1492         __isl_take isl_pw_aff *pwaff2)
1493 {
1494         return pw_aff_union_opt(pwaff1, pwaff2, &isl_aff_ge_basic_set);
1495 }
1496
1497 __isl_give isl_pw_aff *isl_pw_aff_union_max(__isl_take isl_pw_aff *pwaff1,
1498         __isl_take isl_pw_aff *pwaff2)
1499 {
1500         return isl_pw_aff_align_params_pw_pw_and(pwaff1, pwaff2,
1501                                                         &pw_aff_union_max);
1502 }
1503
1504 /* Compute a piecewise quasi-affine expression with a domain that
1505  * is the union of those of pwaff1 and pwaff2 and such that on each
1506  * cell, the quasi-affine expression is the minimum of those of pwaff1
1507  * and pwaff2.  If only one of pwaff1 or pwaff2 is defined on a given
1508  * cell, then the associated expression is the defined one.
1509  */
1510 static __isl_give isl_pw_aff *pw_aff_union_min(__isl_take isl_pw_aff *pwaff1,
1511         __isl_take isl_pw_aff *pwaff2)
1512 {
1513         return pw_aff_union_opt(pwaff1, pwaff2, &isl_aff_le_basic_set);
1514 }
1515
1516 __isl_give isl_pw_aff *isl_pw_aff_union_min(__isl_take isl_pw_aff *pwaff1,
1517         __isl_take isl_pw_aff *pwaff2)
1518 {
1519         return isl_pw_aff_align_params_pw_pw_and(pwaff1, pwaff2,
1520                                                         &pw_aff_union_min);
1521 }
1522
1523 __isl_give isl_pw_aff *isl_pw_aff_union_opt(__isl_take isl_pw_aff *pwaff1,
1524         __isl_take isl_pw_aff *pwaff2, int max)
1525 {
1526         if (max)
1527                 return isl_pw_aff_union_max(pwaff1, pwaff2);
1528         else
1529                 return isl_pw_aff_union_min(pwaff1, pwaff2);
1530 }
1531
1532 /* Construct a map with as domain the domain of pwaff and
1533  * one-dimensional range corresponding to the affine expressions.
1534  */
1535 static __isl_give isl_map *map_from_pw_aff(__isl_take isl_pw_aff *pwaff)
1536 {
1537         int i;
1538         isl_space *dim;
1539         isl_map *map;
1540
1541         if (!pwaff)
1542                 return NULL;
1543
1544         dim = isl_pw_aff_get_space(pwaff);
1545         map = isl_map_empty(dim);
1546
1547         for (i = 0; i < pwaff->n; ++i) {
1548                 isl_basic_map *bmap;
1549                 isl_map *map_i;
1550
1551                 bmap = isl_basic_map_from_aff(isl_aff_copy(pwaff->p[i].aff));
1552                 map_i = isl_map_from_basic_map(bmap);
1553                 map_i = isl_map_intersect_domain(map_i,
1554                                                 isl_set_copy(pwaff->p[i].set));
1555                 map = isl_map_union_disjoint(map, map_i);
1556         }
1557
1558         isl_pw_aff_free(pwaff);
1559
1560         return map;
1561 }
1562
1563 /* Construct a map with as domain the domain of pwaff and
1564  * one-dimensional range corresponding to the affine expressions.
1565  */
1566 __isl_give isl_map *isl_map_from_pw_aff(__isl_take isl_pw_aff *pwaff)
1567 {
1568         if (!pwaff)
1569                 return NULL;
1570         if (isl_space_is_set(pwaff->dim))
1571                 isl_die(isl_pw_aff_get_ctx(pwaff), isl_error_invalid,
1572                         "space of input is not a map",
1573                         return isl_pw_aff_free(pwaff));
1574         return map_from_pw_aff(pwaff);
1575 }
1576
1577 /* Construct a one-dimensional set with as parameter domain
1578  * the domain of pwaff and the single set dimension
1579  * corresponding to the affine expressions.
1580  */
1581 __isl_give isl_set *isl_set_from_pw_aff(__isl_take isl_pw_aff *pwaff)
1582 {
1583         if (!pwaff)
1584                 return NULL;
1585         if (!isl_space_is_set(pwaff->dim))
1586                 isl_die(isl_pw_aff_get_ctx(pwaff), isl_error_invalid,
1587                         "space of input is not a set",
1588                         return isl_pw_aff_free(pwaff));
1589         return map_from_pw_aff(pwaff);
1590 }
1591
1592 /* Return a set containing those elements in the domain
1593  * of pwaff where it is non-negative.
1594  */
1595 __isl_give isl_set *isl_pw_aff_nonneg_set(__isl_take isl_pw_aff *pwaff)
1596 {
1597         int i;
1598         isl_set *set;
1599
1600         if (!pwaff)
1601                 return NULL;
1602
1603         set = isl_set_empty(isl_pw_aff_get_domain_space(pwaff));
1604
1605         for (i = 0; i < pwaff->n; ++i) {
1606                 isl_basic_set *bset;
1607                 isl_set *set_i;
1608
1609                 bset = isl_aff_nonneg_basic_set(isl_aff_copy(pwaff->p[i].aff));
1610                 set_i = isl_set_from_basic_set(bset);
1611                 set_i = isl_set_intersect(set_i, isl_set_copy(pwaff->p[i].set));
1612                 set = isl_set_union_disjoint(set, set_i);
1613         }
1614
1615         isl_pw_aff_free(pwaff);
1616
1617         return set;
1618 }
1619
1620 /* Return a set containing those elements in the domain
1621  * of pwaff where it is zero (if complement is 0) or not zero
1622  * (if complement is 1).
1623  */
1624 static __isl_give isl_set *pw_aff_zero_set(__isl_take isl_pw_aff *pwaff,
1625         int complement)
1626 {
1627         int i;
1628         isl_set *set;
1629
1630         if (!pwaff)
1631                 return NULL;
1632
1633         set = isl_set_empty(isl_pw_aff_get_domain_space(pwaff));
1634
1635         for (i = 0; i < pwaff->n; ++i) {
1636                 isl_basic_set *bset;
1637                 isl_set *set_i, *zero;
1638
1639                 bset = isl_aff_zero_basic_set(isl_aff_copy(pwaff->p[i].aff));
1640                 zero = isl_set_from_basic_set(bset);
1641                 set_i = isl_set_copy(pwaff->p[i].set);
1642                 if (complement)
1643                         set_i = isl_set_subtract(set_i, zero);
1644                 else
1645                         set_i = isl_set_intersect(set_i, zero);
1646                 set = isl_set_union_disjoint(set, set_i);
1647         }
1648
1649         isl_pw_aff_free(pwaff);
1650
1651         return set;
1652 }
1653
1654 /* Return a set containing those elements in the domain
1655  * of pwaff where it is zero.
1656  */
1657 __isl_give isl_set *isl_pw_aff_zero_set(__isl_take isl_pw_aff *pwaff)
1658 {
1659         return pw_aff_zero_set(pwaff, 0);
1660 }
1661
1662 /* Return a set containing those elements in the domain
1663  * of pwaff where it is not zero.
1664  */
1665 __isl_give isl_set *isl_pw_aff_non_zero_set(__isl_take isl_pw_aff *pwaff)
1666 {
1667         return pw_aff_zero_set(pwaff, 1);
1668 }
1669
1670 /* Return a set containing those elements in the shared domain
1671  * of pwaff1 and pwaff2 where pwaff1 is greater than (or equal) to pwaff2.
1672  *
1673  * We compute the difference on the shared domain and then construct
1674  * the set of values where this difference is non-negative.
1675  * If strict is set, we first subtract 1 from the difference.
1676  * If equal is set, we only return the elements where pwaff1 and pwaff2
1677  * are equal.
1678  */
1679 static __isl_give isl_set *pw_aff_gte_set(__isl_take isl_pw_aff *pwaff1,
1680         __isl_take isl_pw_aff *pwaff2, int strict, int equal)
1681 {
1682         isl_set *set1, *set2;
1683
1684         set1 = isl_pw_aff_domain(isl_pw_aff_copy(pwaff1));
1685         set2 = isl_pw_aff_domain(isl_pw_aff_copy(pwaff2));
1686         set1 = isl_set_intersect(set1, set2);
1687         pwaff1 = isl_pw_aff_intersect_domain(pwaff1, isl_set_copy(set1));
1688         pwaff2 = isl_pw_aff_intersect_domain(pwaff2, isl_set_copy(set1));
1689         pwaff1 = isl_pw_aff_add(pwaff1, isl_pw_aff_neg(pwaff2));
1690
1691         if (strict) {
1692                 isl_space *dim = isl_set_get_space(set1);
1693                 isl_aff *aff;
1694                 aff = isl_aff_zero_on_domain(isl_local_space_from_space(dim));
1695                 aff = isl_aff_add_constant_si(aff, -1);
1696                 pwaff1 = isl_pw_aff_add(pwaff1, isl_pw_aff_alloc(set1, aff));
1697         } else
1698                 isl_set_free(set1);
1699
1700         if (equal)
1701                 return isl_pw_aff_zero_set(pwaff1);
1702         return isl_pw_aff_nonneg_set(pwaff1);
1703 }
1704
1705 /* Return a set containing those elements in the shared domain
1706  * of pwaff1 and pwaff2 where pwaff1 is equal to pwaff2.
1707  */
1708 static __isl_give isl_set *pw_aff_eq_set(__isl_take isl_pw_aff *pwaff1,
1709         __isl_take isl_pw_aff *pwaff2)
1710 {
1711         return pw_aff_gte_set(pwaff1, pwaff2, 0, 1);
1712 }
1713
1714 __isl_give isl_set *isl_pw_aff_eq_set(__isl_take isl_pw_aff *pwaff1,
1715         __isl_take isl_pw_aff *pwaff2)
1716 {
1717         return align_params_pw_pw_set_and(pwaff1, pwaff2, &pw_aff_eq_set);
1718 }
1719
1720 /* Return a set containing those elements in the shared domain
1721  * of pwaff1 and pwaff2 where pwaff1 is greater than or equal to pwaff2.
1722  */
1723 static __isl_give isl_set *pw_aff_ge_set(__isl_take isl_pw_aff *pwaff1,
1724         __isl_take isl_pw_aff *pwaff2)
1725 {
1726         return pw_aff_gte_set(pwaff1, pwaff2, 0, 0);
1727 }
1728
1729 __isl_give isl_set *isl_pw_aff_ge_set(__isl_take isl_pw_aff *pwaff1,
1730         __isl_take isl_pw_aff *pwaff2)
1731 {
1732         return align_params_pw_pw_set_and(pwaff1, pwaff2, &pw_aff_ge_set);
1733 }
1734
1735 /* Return a set containing those elements in the shared domain
1736  * of pwaff1 and pwaff2 where pwaff1 is strictly greater than pwaff2.
1737  */
1738 static __isl_give isl_set *pw_aff_gt_set(__isl_take isl_pw_aff *pwaff1,
1739         __isl_take isl_pw_aff *pwaff2)
1740 {
1741         return pw_aff_gte_set(pwaff1, pwaff2, 1, 0);
1742 }
1743
1744 __isl_give isl_set *isl_pw_aff_gt_set(__isl_take isl_pw_aff *pwaff1,
1745         __isl_take isl_pw_aff *pwaff2)
1746 {
1747         return align_params_pw_pw_set_and(pwaff1, pwaff2, &pw_aff_gt_set);
1748 }
1749
1750 __isl_give isl_set *isl_pw_aff_le_set(__isl_take isl_pw_aff *pwaff1,
1751         __isl_take isl_pw_aff *pwaff2)
1752 {
1753         return isl_pw_aff_ge_set(pwaff2, pwaff1);
1754 }
1755
1756 __isl_give isl_set *isl_pw_aff_lt_set(__isl_take isl_pw_aff *pwaff1,
1757         __isl_take isl_pw_aff *pwaff2)
1758 {
1759         return isl_pw_aff_gt_set(pwaff2, pwaff1);
1760 }
1761
1762 /* Return a set containing those elements in the shared domain
1763  * of the elements of list1 and list2 where each element in list1
1764  * has the relation specified by "fn" with each element in list2.
1765  */
1766 static __isl_give isl_set *pw_aff_list_set(__isl_take isl_pw_aff_list *list1,
1767         __isl_take isl_pw_aff_list *list2,
1768         __isl_give isl_set *(*fn)(__isl_take isl_pw_aff *pwaff1,
1769                                     __isl_take isl_pw_aff *pwaff2))
1770 {
1771         int i, j;
1772         isl_ctx *ctx;
1773         isl_set *set;
1774
1775         if (!list1 || !list2)
1776                 goto error;
1777
1778         ctx = isl_pw_aff_list_get_ctx(list1);
1779         if (list1->n < 1 || list2->n < 1)
1780                 isl_die(ctx, isl_error_invalid,
1781                         "list should contain at least one element", goto error);
1782
1783         set = isl_set_universe(isl_pw_aff_get_domain_space(list1->p[0]));
1784         for (i = 0; i < list1->n; ++i)
1785                 for (j = 0; j < list2->n; ++j) {
1786                         isl_set *set_ij;
1787
1788                         set_ij = fn(isl_pw_aff_copy(list1->p[i]),
1789                                     isl_pw_aff_copy(list2->p[j]));
1790                         set = isl_set_intersect(set, set_ij);
1791                 }
1792
1793         isl_pw_aff_list_free(list1);
1794         isl_pw_aff_list_free(list2);
1795         return set;
1796 error:
1797         isl_pw_aff_list_free(list1);
1798         isl_pw_aff_list_free(list2);
1799         return NULL;
1800 }
1801
1802 /* Return a set containing those elements in the shared domain
1803  * of the elements of list1 and list2 where each element in list1
1804  * is equal to each element in list2.
1805  */
1806 __isl_give isl_set *isl_pw_aff_list_eq_set(__isl_take isl_pw_aff_list *list1,
1807         __isl_take isl_pw_aff_list *list2)
1808 {
1809         return pw_aff_list_set(list1, list2, &isl_pw_aff_eq_set);
1810 }
1811
1812 __isl_give isl_set *isl_pw_aff_list_ne_set(__isl_take isl_pw_aff_list *list1,
1813         __isl_take isl_pw_aff_list *list2)
1814 {
1815         return pw_aff_list_set(list1, list2, &isl_pw_aff_ne_set);
1816 }
1817
1818 /* Return a set containing those elements in the shared domain
1819  * of the elements of list1 and list2 where each element in list1
1820  * is less than or equal to each element in list2.
1821  */
1822 __isl_give isl_set *isl_pw_aff_list_le_set(__isl_take isl_pw_aff_list *list1,
1823         __isl_take isl_pw_aff_list *list2)
1824 {
1825         return pw_aff_list_set(list1, list2, &isl_pw_aff_le_set);
1826 }
1827
1828 __isl_give isl_set *isl_pw_aff_list_lt_set(__isl_take isl_pw_aff_list *list1,
1829         __isl_take isl_pw_aff_list *list2)
1830 {
1831         return pw_aff_list_set(list1, list2, &isl_pw_aff_lt_set);
1832 }
1833
1834 __isl_give isl_set *isl_pw_aff_list_ge_set(__isl_take isl_pw_aff_list *list1,
1835         __isl_take isl_pw_aff_list *list2)
1836 {
1837         return pw_aff_list_set(list1, list2, &isl_pw_aff_ge_set);
1838 }
1839
1840 __isl_give isl_set *isl_pw_aff_list_gt_set(__isl_take isl_pw_aff_list *list1,
1841         __isl_take isl_pw_aff_list *list2)
1842 {
1843         return pw_aff_list_set(list1, list2, &isl_pw_aff_gt_set);
1844 }
1845
1846
1847 /* Return a set containing those elements in the shared domain
1848  * of pwaff1 and pwaff2 where pwaff1 is not equal to pwaff2.
1849  */
1850 static __isl_give isl_set *pw_aff_ne_set(__isl_take isl_pw_aff *pwaff1,
1851         __isl_take isl_pw_aff *pwaff2)
1852 {
1853         isl_set *set_lt, *set_gt;
1854
1855         set_lt = isl_pw_aff_lt_set(isl_pw_aff_copy(pwaff1),
1856                                    isl_pw_aff_copy(pwaff2));
1857         set_gt = isl_pw_aff_gt_set(pwaff1, pwaff2);
1858         return isl_set_union_disjoint(set_lt, set_gt);
1859 }
1860
1861 __isl_give isl_set *isl_pw_aff_ne_set(__isl_take isl_pw_aff *pwaff1,
1862         __isl_take isl_pw_aff *pwaff2)
1863 {
1864         return align_params_pw_pw_set_and(pwaff1, pwaff2, &pw_aff_ne_set);
1865 }
1866
1867 __isl_give isl_pw_aff *isl_pw_aff_scale_down(__isl_take isl_pw_aff *pwaff,
1868         isl_int v)
1869 {
1870         int i;
1871
1872         if (isl_int_is_one(v))
1873                 return pwaff;
1874         if (!isl_int_is_pos(v))
1875                 isl_die(isl_pw_aff_get_ctx(pwaff), isl_error_invalid,
1876                         "factor needs to be positive",
1877                         return isl_pw_aff_free(pwaff));
1878         pwaff = isl_pw_aff_cow(pwaff);
1879         if (!pwaff)
1880                 return NULL;
1881         if (pwaff->n == 0)
1882                 return pwaff;
1883
1884         for (i = 0; i < pwaff->n; ++i) {
1885                 pwaff->p[i].aff = isl_aff_scale_down(pwaff->p[i].aff, v);
1886                 if (!pwaff->p[i].aff)
1887                         return isl_pw_aff_free(pwaff);
1888         }
1889
1890         return pwaff;
1891 }
1892
1893 __isl_give isl_pw_aff *isl_pw_aff_floor(__isl_take isl_pw_aff *pwaff)
1894 {
1895         int i;
1896
1897         pwaff = isl_pw_aff_cow(pwaff);
1898         if (!pwaff)
1899                 return NULL;
1900         if (pwaff->n == 0)
1901                 return pwaff;
1902
1903         for (i = 0; i < pwaff->n; ++i) {
1904                 pwaff->p[i].aff = isl_aff_floor(pwaff->p[i].aff);
1905                 if (!pwaff->p[i].aff)
1906                         return isl_pw_aff_free(pwaff);
1907         }
1908
1909         return pwaff;
1910 }
1911
1912 __isl_give isl_pw_aff *isl_pw_aff_ceil(__isl_take isl_pw_aff *pwaff)
1913 {
1914         int i;
1915
1916         pwaff = isl_pw_aff_cow(pwaff);
1917         if (!pwaff)
1918                 return NULL;
1919         if (pwaff->n == 0)
1920                 return pwaff;
1921
1922         for (i = 0; i < pwaff->n; ++i) {
1923                 pwaff->p[i].aff = isl_aff_ceil(pwaff->p[i].aff);
1924                 if (!pwaff->p[i].aff)
1925                         return isl_pw_aff_free(pwaff);
1926         }
1927
1928         return pwaff;
1929 }
1930
1931 /* Assuming that "cond1" and "cond2" are disjoint,
1932  * return an affine expression that is equal to pwaff1 on cond1
1933  * and to pwaff2 on cond2.
1934  */
1935 static __isl_give isl_pw_aff *isl_pw_aff_select(
1936         __isl_take isl_set *cond1, __isl_take isl_pw_aff *pwaff1,
1937         __isl_take isl_set *cond2, __isl_take isl_pw_aff *pwaff2)
1938 {
1939         pwaff1 = isl_pw_aff_intersect_domain(pwaff1, cond1);
1940         pwaff2 = isl_pw_aff_intersect_domain(pwaff2, cond2);
1941
1942         return isl_pw_aff_add_disjoint(pwaff1, pwaff2);
1943 }
1944
1945 /* Return an affine expression that is equal to pwaff_true for elements
1946  * where "cond" is non-zero and to pwaff_false for elements where "cond"
1947  * is zero.
1948  * That is, return cond ? pwaff_true : pwaff_false;
1949  */
1950 __isl_give isl_pw_aff *isl_pw_aff_cond(__isl_take isl_pw_aff *cond,
1951         __isl_take isl_pw_aff *pwaff_true, __isl_take isl_pw_aff *pwaff_false)
1952 {
1953         isl_set *cond_true, *cond_false;
1954
1955         cond_true = isl_pw_aff_non_zero_set(isl_pw_aff_copy(cond));
1956         cond_false = isl_pw_aff_zero_set(cond);
1957         return isl_pw_aff_select(cond_true, pwaff_true,
1958                                  cond_false, pwaff_false);
1959 }
1960
1961 int isl_aff_is_cst(__isl_keep isl_aff *aff)
1962 {
1963         if (!aff)
1964                 return -1;
1965
1966         return isl_seq_first_non_zero(aff->v->el + 2, aff->v->size - 2) == -1;
1967 }
1968
1969 /* Check whether pwaff is a piecewise constant.
1970  */
1971 int isl_pw_aff_is_cst(__isl_keep isl_pw_aff *pwaff)
1972 {
1973         int i;
1974
1975         if (!pwaff)
1976                 return -1;
1977
1978         for (i = 0; i < pwaff->n; ++i) {
1979                 int is_cst = isl_aff_is_cst(pwaff->p[i].aff);
1980                 if (is_cst < 0 || !is_cst)
1981                         return is_cst;
1982         }
1983
1984         return 1;
1985 }
1986
1987 __isl_give isl_aff *isl_aff_mul(__isl_take isl_aff *aff1,
1988         __isl_take isl_aff *aff2)
1989 {
1990         if (!isl_aff_is_cst(aff2) && isl_aff_is_cst(aff1))
1991                 return isl_aff_mul(aff2, aff1);
1992
1993         if (!isl_aff_is_cst(aff2))
1994                 isl_die(isl_aff_get_ctx(aff1), isl_error_invalid,
1995                         "at least one affine expression should be constant",
1996                         goto error);
1997
1998         aff1 = isl_aff_cow(aff1);
1999         if (!aff1 || !aff2)
2000                 goto error;
2001
2002         aff1 = isl_aff_scale(aff1, aff2->v->el[1]);
2003         aff1 = isl_aff_scale_down(aff1, aff2->v->el[0]);
2004
2005         isl_aff_free(aff2);
2006         return aff1;
2007 error:
2008         isl_aff_free(aff1);
2009         isl_aff_free(aff2);
2010         return NULL;
2011 }
2012
2013 static __isl_give isl_pw_aff *pw_aff_add(__isl_take isl_pw_aff *pwaff1,
2014         __isl_take isl_pw_aff *pwaff2)
2015 {
2016         return isl_pw_aff_on_shared_domain(pwaff1, pwaff2, &isl_aff_add);
2017 }
2018
2019 __isl_give isl_pw_aff *isl_pw_aff_add(__isl_take isl_pw_aff *pwaff1,
2020         __isl_take isl_pw_aff *pwaff2)
2021 {
2022         return isl_pw_aff_align_params_pw_pw_and(pwaff1, pwaff2, &pw_aff_add);
2023 }
2024
2025 __isl_give isl_pw_aff *isl_pw_aff_union_add(__isl_take isl_pw_aff *pwaff1,
2026         __isl_take isl_pw_aff *pwaff2)
2027 {
2028         return isl_pw_aff_union_add_(pwaff1, pwaff2);
2029 }
2030
2031 static __isl_give isl_pw_aff *pw_aff_mul(__isl_take isl_pw_aff *pwaff1,
2032         __isl_take isl_pw_aff *pwaff2)
2033 {
2034         return isl_pw_aff_on_shared_domain(pwaff1, pwaff2, &isl_aff_mul);
2035 }
2036
2037 __isl_give isl_pw_aff *isl_pw_aff_mul(__isl_take isl_pw_aff *pwaff1,
2038         __isl_take isl_pw_aff *pwaff2)
2039 {
2040         return isl_pw_aff_align_params_pw_pw_and(pwaff1, pwaff2, &pw_aff_mul);
2041 }
2042
2043 static __isl_give isl_pw_aff *pw_aff_min(__isl_take isl_pw_aff *pwaff1,
2044         __isl_take isl_pw_aff *pwaff2)
2045 {
2046         isl_set *le;
2047         isl_set *dom;
2048
2049         dom = isl_set_intersect(isl_pw_aff_domain(isl_pw_aff_copy(pwaff1)),
2050                                 isl_pw_aff_domain(isl_pw_aff_copy(pwaff2)));
2051         le = isl_pw_aff_le_set(isl_pw_aff_copy(pwaff1),
2052                                 isl_pw_aff_copy(pwaff2));
2053         dom = isl_set_subtract(dom, isl_set_copy(le));
2054         return isl_pw_aff_select(le, pwaff1, dom, pwaff2);
2055 }
2056
2057 __isl_give isl_pw_aff *isl_pw_aff_min(__isl_take isl_pw_aff *pwaff1,
2058         __isl_take isl_pw_aff *pwaff2)
2059 {
2060         return isl_pw_aff_align_params_pw_pw_and(pwaff1, pwaff2, &pw_aff_min);
2061 }
2062
2063 static __isl_give isl_pw_aff *pw_aff_max(__isl_take isl_pw_aff *pwaff1,
2064         __isl_take isl_pw_aff *pwaff2)
2065 {
2066         isl_set *ge;
2067         isl_set *dom;
2068
2069         dom = isl_set_intersect(isl_pw_aff_domain(isl_pw_aff_copy(pwaff1)),
2070                                 isl_pw_aff_domain(isl_pw_aff_copy(pwaff2)));
2071         ge = isl_pw_aff_ge_set(isl_pw_aff_copy(pwaff1),
2072                                 isl_pw_aff_copy(pwaff2));
2073         dom = isl_set_subtract(dom, isl_set_copy(ge));
2074         return isl_pw_aff_select(ge, pwaff1, dom, pwaff2);
2075 }
2076
2077 __isl_give isl_pw_aff *isl_pw_aff_max(__isl_take isl_pw_aff *pwaff1,
2078         __isl_take isl_pw_aff *pwaff2)
2079 {
2080         return isl_pw_aff_align_params_pw_pw_and(pwaff1, pwaff2, &pw_aff_max);
2081 }
2082
2083 static __isl_give isl_pw_aff *pw_aff_list_reduce(
2084         __isl_take isl_pw_aff_list *list,
2085         __isl_give isl_pw_aff *(*fn)(__isl_take isl_pw_aff *pwaff1,
2086                                         __isl_take isl_pw_aff *pwaff2))
2087 {
2088         int i;
2089         isl_ctx *ctx;
2090         isl_pw_aff *res;
2091
2092         if (!list)
2093                 return NULL;
2094
2095         ctx = isl_pw_aff_list_get_ctx(list);
2096         if (list->n < 1)
2097                 isl_die(ctx, isl_error_invalid,
2098                         "list should contain at least one element",
2099                         return isl_pw_aff_list_free(list));
2100
2101         res = isl_pw_aff_copy(list->p[0]);
2102         for (i = 1; i < list->n; ++i)
2103                 res = fn(res, isl_pw_aff_copy(list->p[i]));
2104
2105         isl_pw_aff_list_free(list);
2106         return res;
2107 }
2108
2109 /* Return an isl_pw_aff that maps each element in the intersection of the
2110  * domains of the elements of list to the minimal corresponding affine
2111  * expression.
2112  */
2113 __isl_give isl_pw_aff *isl_pw_aff_list_min(__isl_take isl_pw_aff_list *list)
2114 {
2115         return pw_aff_list_reduce(list, &isl_pw_aff_min);
2116 }
2117
2118 /* Return an isl_pw_aff that maps each element in the intersection of the
2119  * domains of the elements of list to the maximal corresponding affine
2120  * expression.
2121  */
2122 __isl_give isl_pw_aff *isl_pw_aff_list_max(__isl_take isl_pw_aff_list *list)
2123 {
2124         return pw_aff_list_reduce(list, &isl_pw_aff_max);
2125 }
2126
2127 #undef BASE
2128 #define BASE aff
2129
2130 #include <isl_multi_templ.c>
2131
2132 /* Construct an isl_multi_aff in the given space with value zero in
2133  * each of the output dimensions.
2134  */
2135 __isl_give isl_multi_aff *isl_multi_aff_zero(__isl_take isl_space *space)
2136 {
2137         int n;
2138         isl_multi_aff *ma;
2139
2140         if (!space)
2141                 return NULL;
2142
2143         n = isl_space_dim(space , isl_dim_out);
2144         ma = isl_multi_aff_alloc(isl_space_copy(space));
2145
2146         if (!n)
2147                 isl_space_free(space);
2148         else {
2149                 int i;
2150                 isl_local_space *ls;
2151                 isl_aff *aff;
2152
2153                 space = isl_space_domain(space);
2154                 ls = isl_local_space_from_space(space);
2155                 aff = isl_aff_zero_on_domain(ls);
2156
2157                 for (i = 0; i < n; ++i)
2158                         ma = isl_multi_aff_set_aff(ma, i, isl_aff_copy(aff));
2159
2160                 isl_aff_free(aff);
2161         }
2162
2163         return ma;
2164 }
2165
2166 /* Create an isl_pw_multi_aff with the given isl_multi_aff on a universe
2167  * domain.
2168  */
2169 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_from_multi_aff(
2170         __isl_take isl_multi_aff *ma)
2171 {
2172         isl_set *dom = isl_set_universe(isl_multi_aff_get_domain_space(ma));
2173         return isl_pw_multi_aff_alloc(dom, ma);
2174 }
2175
2176 __isl_give isl_multi_aff *isl_multi_aff_add(__isl_take isl_multi_aff *maff1,
2177         __isl_take isl_multi_aff *maff2)
2178 {
2179         int i;
2180         isl_ctx *ctx;
2181
2182         maff1 = isl_multi_aff_cow(maff1);
2183         if (!maff1 || !maff2)
2184                 goto error;
2185
2186         ctx = isl_multi_aff_get_ctx(maff1);
2187         if (!isl_space_is_equal(maff1->space, maff2->space))
2188                 isl_die(ctx, isl_error_invalid,
2189                         "spaces don't match", goto error);
2190
2191         for (i = 0; i < maff1->n; ++i) {
2192                 maff1->p[i] = isl_aff_add(maff1->p[i],
2193                                             isl_aff_copy(maff2->p[i]));
2194                 if (!maff1->p[i])
2195                         goto error;
2196         }
2197
2198         isl_multi_aff_free(maff2);
2199         return maff1;
2200 error:
2201         isl_multi_aff_free(maff1);
2202         isl_multi_aff_free(maff2);
2203         return NULL;
2204 }
2205
2206 /* Exploit the equalities in "eq" to simplify the affine expressions.
2207  */
2208 static __isl_give isl_multi_aff *isl_multi_aff_substitute_equalities(
2209         __isl_take isl_multi_aff *maff, __isl_take isl_basic_set *eq)
2210 {
2211         int i;
2212
2213         maff = isl_multi_aff_cow(maff);
2214         if (!maff || !eq)
2215                 goto error;
2216
2217         for (i = 0; i < maff->n; ++i) {
2218                 maff->p[i] = isl_aff_substitute_equalities(maff->p[i],
2219                                                     isl_basic_set_copy(eq));
2220                 if (!maff->p[i])
2221                         goto error;
2222         }
2223
2224         isl_basic_set_free(eq);
2225         return maff;
2226 error:
2227         isl_basic_set_free(eq);
2228         isl_multi_aff_free(maff);
2229         return NULL;
2230 }
2231
2232 __isl_give isl_multi_aff *isl_multi_aff_scale(__isl_take isl_multi_aff *maff,
2233         isl_int f)
2234 {
2235         int i;
2236
2237         maff = isl_multi_aff_cow(maff);
2238         if (!maff)
2239                 return NULL;
2240
2241         for (i = 0; i < maff->n; ++i) {
2242                 maff->p[i] = isl_aff_scale(maff->p[i], f);
2243                 if (!maff->p[i])
2244                         return isl_multi_aff_free(maff);
2245         }
2246
2247         return maff;
2248 }
2249
2250 __isl_give isl_multi_aff *isl_multi_aff_add_on_domain(__isl_keep isl_set *dom,
2251         __isl_take isl_multi_aff *maff1, __isl_take isl_multi_aff *maff2)
2252 {
2253         maff1 = isl_multi_aff_add(maff1, maff2);
2254         maff1 = isl_multi_aff_gist(maff1, isl_set_copy(dom));
2255         return maff1;
2256 }
2257
2258 int isl_multi_aff_is_empty(__isl_keep isl_multi_aff *maff)
2259 {
2260         if (!maff)
2261                 return -1;
2262
2263         return 0;
2264 }
2265
2266 int isl_multi_aff_plain_is_equal(__isl_keep isl_multi_aff *maff1,
2267         __isl_keep isl_multi_aff *maff2)
2268 {
2269         int i;
2270         int equal;
2271
2272         if (!maff1 || !maff2)
2273                 return -1;
2274         if (maff1->n != maff2->n)
2275                 return 0;
2276         equal = isl_space_is_equal(maff1->space, maff2->space);
2277         if (equal < 0 || !equal)
2278                 return equal;
2279
2280         for (i = 0; i < maff1->n; ++i) {
2281                 equal = isl_aff_plain_is_equal(maff1->p[i], maff2->p[i]);
2282                 if (equal < 0 || !equal)
2283                         return equal;
2284         }
2285
2286         return 1;
2287 }
2288
2289 __isl_give isl_multi_aff *isl_multi_aff_set_dim_name(
2290         __isl_take isl_multi_aff *maff,
2291         enum isl_dim_type type, unsigned pos, const char *s)
2292 {
2293         int i;
2294
2295         maff = isl_multi_aff_cow(maff);
2296         if (!maff)
2297                 return NULL;
2298
2299         maff->space = isl_space_set_dim_name(maff->space, type, pos, s);
2300         if (!maff->space)
2301                 return isl_multi_aff_free(maff);
2302
2303         if (type == isl_dim_out)
2304                 return maff;
2305         for (i = 0; i < maff->n; ++i) {
2306                 maff->p[i] = isl_aff_set_dim_name(maff->p[i], type, pos, s);
2307                 if (!maff->p[i])
2308                         return isl_multi_aff_free(maff);
2309         }
2310
2311         return maff;
2312 }
2313
2314 __isl_give isl_multi_aff *isl_multi_aff_drop_dims(__isl_take isl_multi_aff *maff,
2315         enum isl_dim_type type, unsigned first, unsigned n)
2316 {
2317         int i;
2318
2319         maff = isl_multi_aff_cow(maff);
2320         if (!maff)
2321                 return NULL;
2322
2323         maff->space = isl_space_drop_dims(maff->space, type, first, n);
2324         if (!maff->space)
2325                 return isl_multi_aff_free(maff);
2326
2327         if (type == isl_dim_out) {
2328                 for (i = 0; i < n; ++i)
2329                         isl_aff_free(maff->p[first + i]);
2330                 for (i = first; i + n < maff->n; ++i)
2331                         maff->p[i] = maff->p[i + n];
2332                 maff->n -= n;
2333                 return maff;
2334         }
2335
2336         for (i = 0; i < maff->n; ++i) {
2337                 maff->p[i] = isl_aff_drop_dims(maff->p[i], type, first, n);
2338                 if (!maff->p[i])
2339                         return isl_multi_aff_free(maff);
2340         }
2341
2342         return maff;
2343 }
2344
2345 #undef PW
2346 #define PW isl_pw_multi_aff
2347 #undef EL
2348 #define EL isl_multi_aff
2349 #undef EL_IS_ZERO
2350 #define EL_IS_ZERO is_empty
2351 #undef ZERO
2352 #define ZERO empty
2353 #undef IS_ZERO
2354 #define IS_ZERO is_empty
2355 #undef FIELD
2356 #define FIELD maff
2357 #undef DEFAULT_IS_ZERO
2358 #define DEFAULT_IS_ZERO 0
2359
2360 #define NO_NEG
2361 #define NO_EVAL
2362 #define NO_OPT
2363 #define NO_INVOLVES_DIMS
2364 #define NO_MOVE_DIMS
2365 #define NO_INSERT_DIMS
2366 #define NO_LIFT
2367 #define NO_MORPH
2368
2369 #include <isl_pw_templ.c>
2370
2371 #undef UNION
2372 #define UNION isl_union_pw_multi_aff
2373 #undef PART
2374 #define PART isl_pw_multi_aff
2375 #undef PARTS
2376 #define PARTS pw_multi_aff
2377 #define ALIGN_DOMAIN
2378
2379 #define NO_EVAL
2380
2381 #include <isl_union_templ.c>
2382
2383 static __isl_give isl_pw_multi_aff *pw_multi_aff_add(
2384         __isl_take isl_pw_multi_aff *pma1, __isl_take isl_pw_multi_aff *pma2)
2385 {
2386         return isl_pw_multi_aff_on_shared_domain(pma1, pma2,
2387                                                 &isl_multi_aff_add);
2388 }
2389
2390 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_add(
2391         __isl_take isl_pw_multi_aff *pma1, __isl_take isl_pw_multi_aff *pma2)
2392 {
2393         return isl_pw_multi_aff_align_params_pw_pw_and(pma1, pma2,
2394                                                 &pw_multi_aff_add);
2395 }
2396
2397 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_union_add(
2398         __isl_take isl_pw_multi_aff *pma1, __isl_take isl_pw_multi_aff *pma2)
2399 {
2400         return isl_pw_multi_aff_union_add_(pma1, pma2);
2401 }
2402
2403 /* Construct a map mapping the domain of the piecewise multi-affine expression
2404  * to its range, with each dimension in the range equated to the
2405  * corresponding affine expression on its cell.
2406  */
2407 __isl_give isl_map *isl_map_from_pw_multi_aff(__isl_take isl_pw_multi_aff *pma)
2408 {
2409         int i;
2410         isl_map *map;
2411
2412         if (!pma)
2413                 return NULL;
2414
2415         map = isl_map_empty(isl_pw_multi_aff_get_space(pma));
2416
2417         for (i = 0; i < pma->n; ++i) {
2418                 isl_multi_aff *maff;
2419                 isl_basic_map *bmap;
2420                 isl_map *map_i;
2421
2422                 maff = isl_multi_aff_copy(pma->p[i].maff);
2423                 bmap = isl_basic_map_from_multi_aff(maff);
2424                 map_i = isl_map_from_basic_map(bmap);
2425                 map_i = isl_map_intersect_domain(map_i,
2426                                                 isl_set_copy(pma->p[i].set));
2427                 map = isl_map_union_disjoint(map, map_i);
2428         }
2429
2430         isl_pw_multi_aff_free(pma);
2431         return map;
2432 }
2433
2434 __isl_give isl_set *isl_set_from_pw_multi_aff(__isl_take isl_pw_multi_aff *pma)
2435 {
2436         if (!isl_space_is_set(pma->dim))
2437                 isl_die(isl_pw_multi_aff_get_ctx(pma), isl_error_invalid,
2438                         "isl_pw_multi_aff cannot be converted into an isl_set",
2439                         return isl_pw_multi_aff_free(pma));
2440
2441         return isl_map_from_pw_multi_aff(pma);
2442 }
2443
2444 /* Given a basic map with a single output dimension that is defined
2445  * in terms of the parameters and input dimensions using an equality,
2446  * extract an isl_aff that expresses the output dimension in terms
2447  * of the parameters and input dimensions.
2448  *
2449  * Since some applications expect the result of isl_pw_multi_aff_from_map
2450  * to only contain integer affine expressions, we compute the floor
2451  * of the expression before returning.
2452  *
2453  * This function shares some similarities with
2454  * isl_basic_map_has_defining_equality and isl_constraint_get_bound.
2455  */
2456 static __isl_give isl_aff *extract_isl_aff_from_basic_map(
2457         __isl_take isl_basic_map *bmap)
2458 {
2459         int i;
2460         unsigned offset;
2461         unsigned total;
2462         isl_local_space *ls;
2463         isl_aff *aff;
2464
2465         if (!bmap)
2466                 return NULL;
2467         if (isl_basic_map_dim(bmap, isl_dim_out) != 1)
2468                 isl_die(isl_basic_map_get_ctx(bmap), isl_error_invalid,
2469                         "basic map should have a single output dimension",
2470                         goto error);
2471         offset = isl_basic_map_offset(bmap, isl_dim_out);
2472         total = isl_basic_map_total_dim(bmap);
2473         for (i = 0; i < bmap->n_eq; ++i) {
2474                 if (isl_int_is_zero(bmap->eq[i][offset]))
2475                         continue;
2476                 if (isl_seq_first_non_zero(bmap->eq[i] + offset + 1,
2477                                            1 + total - (offset + 1)) != -1)
2478                         continue;
2479                 break;
2480         }
2481         if (i >= bmap->n_eq)
2482                 isl_die(isl_basic_map_get_ctx(bmap), isl_error_invalid,
2483                         "unable to find suitable equality", goto error);
2484         ls = isl_basic_map_get_local_space(bmap);
2485         aff = isl_aff_alloc(isl_local_space_domain(ls));
2486         if (!aff)
2487                 goto error;
2488         if (isl_int_is_neg(bmap->eq[i][offset]))
2489                 isl_seq_cpy(aff->v->el + 1, bmap->eq[i], offset);
2490         else
2491                 isl_seq_neg(aff->v->el + 1, bmap->eq[i], offset);
2492         isl_seq_clr(aff->v->el + 1 + offset, aff->v->size - (1 + offset));
2493         isl_int_abs(aff->v->el[0], bmap->eq[i][offset]);
2494         isl_basic_map_free(bmap);
2495
2496         aff = isl_aff_remove_unused_divs(aff);
2497         aff = isl_aff_floor(aff);
2498         return aff;
2499 error:
2500         isl_basic_map_free(bmap);
2501         return NULL;
2502 }
2503
2504 /* Given a basic map where each output dimension is defined
2505  * in terms of the parameters and input dimensions using an equality,
2506  * extract an isl_multi_aff that expresses the output dimensions in terms
2507  * of the parameters and input dimensions.
2508  */
2509 static __isl_give isl_multi_aff *extract_isl_multi_aff_from_basic_map(
2510         __isl_take isl_basic_map *bmap)
2511 {
2512         int i;
2513         unsigned n_out;
2514         isl_multi_aff *ma;
2515
2516         if (!bmap)
2517                 return NULL;
2518
2519         ma = isl_multi_aff_alloc(isl_basic_map_get_space(bmap));
2520         n_out = isl_basic_map_dim(bmap, isl_dim_out);
2521
2522         for (i = 0; i < n_out; ++i) {
2523                 isl_basic_map *bmap_i;
2524                 isl_aff *aff;
2525
2526                 bmap_i = isl_basic_map_copy(bmap);
2527                 bmap_i = isl_basic_map_project_out(bmap_i, isl_dim_out,
2528                                                         i + 1, n_out - (1 + i));
2529                 bmap_i = isl_basic_map_project_out(bmap_i, isl_dim_out, 0, i);
2530                 aff = extract_isl_aff_from_basic_map(bmap_i);
2531                 ma = isl_multi_aff_set_aff(ma, i, aff);
2532         }
2533
2534         isl_basic_map_free(bmap);
2535
2536         return ma;
2537 }
2538
2539 /* Create an isl_pw_multi_aff that is equivalent to
2540  * isl_map_intersect_domain(isl_map_from_basic_map(bmap), domain).
2541  * The given basic map is such that each output dimension is defined
2542  * in terms of the parameters and input dimensions using an equality.
2543  */
2544 static __isl_give isl_pw_multi_aff *plain_pw_multi_aff_from_map(
2545         __isl_take isl_set *domain, __isl_take isl_basic_map *bmap)
2546 {
2547         isl_multi_aff *ma;
2548
2549         ma = extract_isl_multi_aff_from_basic_map(bmap);
2550         return isl_pw_multi_aff_alloc(domain, ma);
2551 }
2552
2553 /* Try and create an isl_pw_multi_aff that is equivalent to the given isl_map.
2554  * This obivously only works if the input "map" is single-valued.
2555  * If so, we compute the lexicographic minimum of the image in the form
2556  * of an isl_pw_multi_aff.  Since the image is unique, it is equal
2557  * to its lexicographic minimum.
2558  * If the input is not single-valued, we produce an error.
2559  *
2560  * As a special case, we first check if all output dimensions are uniquely
2561  * defined in terms of the parameters and input dimensions over the entire
2562  * domain.  If so, we extract the desired isl_pw_multi_aff directly
2563  * from the affine hull of "map" and its domain.
2564  */
2565 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_from_map(__isl_take isl_map *map)
2566 {
2567         int i;
2568         int sv;
2569         isl_pw_multi_aff *pma;
2570         isl_basic_map *hull;
2571
2572         if (!map)
2573                 return NULL;
2574
2575         hull = isl_map_affine_hull(isl_map_copy(map));
2576         sv = isl_basic_map_plain_is_single_valued(hull);
2577         if (sv >= 0 && sv)
2578                 return plain_pw_multi_aff_from_map(isl_map_domain(map), hull);
2579         isl_basic_map_free(hull);
2580         if (sv < 0)
2581                 goto error;
2582
2583         sv = isl_map_is_single_valued(map);
2584         if (sv < 0)
2585                 goto error;
2586         if (!sv)
2587                 isl_die(isl_map_get_ctx(map), isl_error_invalid,
2588                         "map is not single-valued", goto error);
2589         map = isl_map_make_disjoint(map);
2590         if (!map)
2591                 return NULL;
2592
2593         pma = isl_pw_multi_aff_empty(isl_map_get_space(map));
2594
2595         for (i = 0; i < map->n; ++i) {
2596                 isl_pw_multi_aff *pma_i;
2597                 isl_basic_map *bmap;
2598                 bmap = isl_basic_map_copy(map->p[i]);
2599                 pma_i = isl_basic_map_lexmin_pw_multi_aff(bmap);
2600                 pma = isl_pw_multi_aff_add_disjoint(pma, pma_i);
2601         }
2602
2603         isl_map_free(map);
2604         return pma;
2605 error:
2606         isl_map_free(map);
2607         return NULL;
2608 }
2609
2610 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_from_set(__isl_take isl_set *set)
2611 {
2612         return isl_pw_multi_aff_from_map(set);
2613 }
2614
2615 /* Return the piecewise affine expression "set ? 1 : 0".
2616  */
2617 __isl_give isl_pw_aff *isl_set_indicator_function(__isl_take isl_set *set)
2618 {
2619         isl_pw_aff *pa;
2620         isl_space *space = isl_set_get_space(set);
2621         isl_local_space *ls = isl_local_space_from_space(space);
2622         isl_aff *zero = isl_aff_zero_on_domain(isl_local_space_copy(ls));
2623         isl_aff *one = isl_aff_zero_on_domain(ls);
2624
2625         one = isl_aff_add_constant_si(one, 1);
2626         pa = isl_pw_aff_alloc(isl_set_copy(set), one);
2627         set = isl_set_complement(set);
2628         pa = isl_pw_aff_add_disjoint(pa, isl_pw_aff_alloc(set, zero));
2629
2630         return pa;
2631 }
2632
2633 /* Plug in "subs" for dimension "type", "pos" of "aff".
2634  *
2635  * Let i be the dimension to replace and let "subs" be of the form
2636  *
2637  *      f/d
2638  *
2639  * and "aff" of the form
2640  *
2641  *      (a i + g)/m
2642  *
2643  * The result is
2644  *
2645  *      floor((a f + d g')/(m d))
2646  *
2647  * where g' is the result of plugging in "subs" in each of the integer
2648  * divisions in g.
2649  */
2650 __isl_give isl_aff *isl_aff_substitute(__isl_take isl_aff *aff,
2651         enum isl_dim_type type, unsigned pos, __isl_keep isl_aff *subs)
2652 {
2653         isl_ctx *ctx;
2654         isl_int v;
2655
2656         aff = isl_aff_cow(aff);
2657         if (!aff || !subs)
2658                 return isl_aff_free(aff);
2659
2660         ctx = isl_aff_get_ctx(aff);
2661         if (!isl_space_is_equal(aff->ls->dim, subs->ls->dim))
2662                 isl_die(ctx, isl_error_invalid,
2663                         "spaces don't match", return isl_aff_free(aff));
2664         if (isl_local_space_dim(subs->ls, isl_dim_div) != 0)
2665                 isl_die(ctx, isl_error_unsupported,
2666                         "cannot handle divs yet", return isl_aff_free(aff));
2667
2668         aff->ls = isl_local_space_substitute(aff->ls, type, pos, subs);
2669         if (!aff->ls)
2670                 return isl_aff_free(aff);
2671
2672         aff->v = isl_vec_cow(aff->v);
2673         if (!aff->v)
2674                 return isl_aff_free(aff);
2675
2676         pos += isl_local_space_offset(aff->ls, type);
2677
2678         isl_int_init(v);
2679         isl_int_set(v, aff->v->el[1 + pos]);
2680         isl_int_set_si(aff->v->el[1 + pos], 0);
2681         isl_seq_combine(aff->v->el + 1, subs->v->el[0], aff->v->el + 1,
2682                         v, subs->v->el + 1, subs->v->size - 1);
2683         isl_int_mul(aff->v->el[0], aff->v->el[0], subs->v->el[0]);
2684         isl_int_clear(v);
2685
2686         return aff;
2687 }
2688
2689 /* Plug in "subs" for dimension "type", "pos" in each of the affine
2690  * expressions in "maff".
2691  */
2692 __isl_give isl_multi_aff *isl_multi_aff_substitute(
2693         __isl_take isl_multi_aff *maff, enum isl_dim_type type, unsigned pos,
2694         __isl_keep isl_aff *subs)
2695 {
2696         int i;
2697
2698         maff = isl_multi_aff_cow(maff);
2699         if (!maff || !subs)
2700                 return isl_multi_aff_free(maff);
2701
2702         if (type == isl_dim_in)
2703                 type = isl_dim_set;
2704
2705         for (i = 0; i < maff->n; ++i) {
2706                 maff->p[i] = isl_aff_substitute(maff->p[i], type, pos, subs);
2707                 if (!maff->p[i])
2708                         return isl_multi_aff_free(maff);
2709         }
2710
2711         return maff;
2712 }
2713
2714 /* Plug in "subs" for dimension "type", "pos" of "pma".
2715  *
2716  * pma is of the form
2717  *
2718  *      A_i(v) -> M_i(v)
2719  *
2720  * while subs is of the form
2721  *
2722  *      v' = B_j(v) -> S_j
2723  *
2724  * Each pair i,j such that C_ij = A_i \cap B_i is non-empty
2725  * has a contribution in the result, in particular
2726  *
2727  *      C_ij(S_j) -> M_i(S_j)
2728  *
2729  * Note that plugging in S_j in C_ij may also result in an empty set
2730  * and this contribution should simply be discarded.
2731  */
2732 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_substitute(
2733         __isl_take isl_pw_multi_aff *pma, enum isl_dim_type type, unsigned pos,
2734         __isl_keep isl_pw_aff *subs)
2735 {
2736         int i, j, n;
2737         isl_pw_multi_aff *res;
2738
2739         if (!pma || !subs)
2740                 return isl_pw_multi_aff_free(pma);
2741
2742         n = pma->n * subs->n;
2743         res = isl_pw_multi_aff_alloc_size(isl_space_copy(pma->dim), n);
2744
2745         for (i = 0; i < pma->n; ++i) {
2746                 for (j = 0; j < subs->n; ++j) {
2747                         isl_set *common;
2748                         isl_multi_aff *res_ij;
2749                         common = isl_set_intersect(
2750                                         isl_set_copy(pma->p[i].set),
2751                                         isl_set_copy(subs->p[j].set));
2752                         common = isl_set_substitute(common,
2753                                         type, pos, subs->p[j].aff);
2754                         if (isl_set_plain_is_empty(common)) {
2755                                 isl_set_free(common);
2756                                 continue;
2757                         }
2758
2759                         res_ij = isl_multi_aff_substitute(
2760                                         isl_multi_aff_copy(pma->p[i].maff),
2761                                         type, pos, subs->p[j].aff);
2762
2763                         res = isl_pw_multi_aff_add_piece(res, common, res_ij);
2764                 }
2765         }
2766
2767         isl_pw_multi_aff_free(pma);
2768         return res;
2769 }
2770
2771 /* Extend the local space of "dst" to include the divs
2772  * in the local space of "src".
2773  */
2774 __isl_give isl_aff *isl_aff_align_divs(__isl_take isl_aff *dst,
2775         __isl_keep isl_aff *src)
2776 {
2777         isl_ctx *ctx;
2778         int *exp1 = NULL;
2779         int *exp2 = NULL;
2780         isl_mat *div;
2781
2782         if (!src || !dst)
2783                 return isl_aff_free(dst);
2784
2785         ctx = isl_aff_get_ctx(src);
2786         if (!isl_space_is_equal(src->ls->dim, dst->ls->dim))
2787                 isl_die(ctx, isl_error_invalid,
2788                         "spaces don't match", goto error);
2789
2790         if (src->ls->div->n_row == 0)
2791                 return dst;
2792
2793         exp1 = isl_alloc_array(ctx, int, src->ls->div->n_row);
2794         exp2 = isl_alloc_array(ctx, int, dst->ls->div->n_row);
2795         if (!exp1 || !exp2)
2796                 goto error;
2797
2798         div = isl_merge_divs(src->ls->div, dst->ls->div, exp1, exp2);
2799         dst = isl_aff_expand_divs(dst, div, exp2);
2800         free(exp1);
2801         free(exp2);
2802
2803         return dst;
2804 error:
2805         free(exp1);
2806         free(exp2);
2807         return isl_aff_free(dst);
2808 }
2809
2810 /* Adjust the local spaces of the affine expressions in "maff"
2811  * such that they all have the save divs.
2812  */
2813 __isl_give isl_multi_aff *isl_multi_aff_align_divs(
2814         __isl_take isl_multi_aff *maff)
2815 {
2816         int i;
2817
2818         if (!maff)
2819                 return NULL;
2820         if (maff->n == 0)
2821                 return maff;
2822         maff = isl_multi_aff_cow(maff);
2823         if (!maff)
2824                 return NULL;
2825
2826         for (i = 1; i < maff->n; ++i)
2827                 maff->p[0] = isl_aff_align_divs(maff->p[0], maff->p[i]);
2828         for (i = 1; i < maff->n; ++i) {
2829                 maff->p[i] = isl_aff_align_divs(maff->p[i], maff->p[0]);
2830                 if (!maff->p[i])
2831                         return isl_multi_aff_free(maff);
2832         }
2833
2834         return maff;
2835 }
2836
2837 __isl_give isl_aff *isl_aff_lift(__isl_take isl_aff *aff)
2838 {
2839         aff = isl_aff_cow(aff);
2840         if (!aff)
2841                 return NULL;
2842
2843         aff->ls = isl_local_space_lift(aff->ls);
2844         if (!aff->ls)
2845                 return isl_aff_free(aff);
2846
2847         return aff;
2848 }
2849
2850 /* Lift "maff" to a space with extra dimensions such that the result
2851  * has no more existentially quantified variables.
2852  * If "ls" is not NULL, then *ls is assigned the local space that lies
2853  * at the basis of the lifting applied to "maff".
2854  */
2855 __isl_give isl_multi_aff *isl_multi_aff_lift(__isl_take isl_multi_aff *maff,
2856         __isl_give isl_local_space **ls)
2857 {
2858         int i;
2859         isl_space *space;
2860         unsigned n_div;
2861
2862         if (ls)
2863                 *ls = NULL;
2864
2865         if (!maff)
2866                 return NULL;
2867
2868         if (maff->n == 0) {
2869                 if (ls) {
2870                         isl_space *space = isl_multi_aff_get_domain_space(maff);
2871                         *ls = isl_local_space_from_space(space);
2872                         if (!*ls)
2873                                 return isl_multi_aff_free(maff);
2874                 }
2875                 return maff;
2876         }
2877
2878         maff = isl_multi_aff_cow(maff);
2879         maff = isl_multi_aff_align_divs(maff);
2880         if (!maff)
2881                 return NULL;
2882
2883         n_div = isl_aff_dim(maff->p[0], isl_dim_div);
2884         space = isl_multi_aff_get_space(maff);
2885         space = isl_space_lift(isl_space_domain(space), n_div);
2886         space = isl_space_extend_domain_with_range(space,
2887                                                 isl_multi_aff_get_space(maff));
2888         if (!space)
2889                 return isl_multi_aff_free(maff);
2890         isl_space_free(maff->space);
2891         maff->space = space;
2892
2893         if (ls) {
2894                 *ls = isl_aff_get_domain_local_space(maff->p[0]);
2895                 if (!*ls)
2896                         return isl_multi_aff_free(maff);
2897         }
2898
2899         for (i = 0; i < maff->n; ++i) {
2900                 maff->p[i] = isl_aff_lift(maff->p[i]);
2901                 if (!maff->p[i])
2902                         goto error;
2903         }
2904
2905         return maff;
2906 error:
2907         if (ls)
2908                 isl_local_space_free(*ls);
2909         return isl_multi_aff_free(maff);
2910 }
2911
2912
2913 /* Extract an isl_pw_aff corresponding to output dimension "pos" of "pma".
2914  */
2915 __isl_give isl_pw_aff *isl_pw_multi_aff_get_pw_aff(
2916         __isl_keep isl_pw_multi_aff *pma, int pos)
2917 {
2918         int i;
2919         int n_out;
2920         isl_space *space;
2921         isl_pw_aff *pa;
2922
2923         if (!pma)
2924                 return NULL;
2925
2926         n_out = isl_pw_multi_aff_dim(pma, isl_dim_out);
2927         if (pos < 0 || pos >= n_out)
2928                 isl_die(isl_pw_multi_aff_get_ctx(pma), isl_error_invalid,
2929                         "index out of bounds", return NULL);
2930
2931         space = isl_pw_multi_aff_get_space(pma);
2932         space = isl_space_drop_dims(space, isl_dim_out,
2933                                     pos + 1, n_out - pos - 1);
2934         space = isl_space_drop_dims(space, isl_dim_out, 0, pos);
2935
2936         pa = isl_pw_aff_alloc_size(space, pma->n);
2937         for (i = 0; i < pma->n; ++i) {
2938                 isl_aff *aff;
2939                 aff = isl_multi_aff_get_aff(pma->p[i].maff, pos);
2940                 pa = isl_pw_aff_add_piece(pa, isl_set_copy(pma->p[i].set), aff);
2941         }
2942
2943         return pa;
2944 }
2945
2946 /* Return an isl_pw_multi_aff with the given "set" as domain and
2947  * an unnamed zero-dimensional range.
2948  */
2949 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_from_domain(
2950         __isl_take isl_set *set)
2951 {
2952         isl_multi_aff *ma;
2953         isl_space *space;
2954
2955         space = isl_set_get_space(set);
2956         space = isl_space_from_domain(space);
2957         ma = isl_multi_aff_zero(space);
2958         return isl_pw_multi_aff_alloc(set, ma);
2959 }
2960
2961 /* Add an isl_pw_multi_aff with the given "set" as domain and
2962  * an unnamed zero-dimensional range to *user.
2963  */
2964 static int add_pw_multi_aff_from_domain(__isl_take isl_set *set, void *user)
2965 {
2966         isl_union_pw_multi_aff **upma = user;
2967         isl_pw_multi_aff *pma;
2968
2969         pma = isl_pw_multi_aff_from_domain(set);
2970         *upma = isl_union_pw_multi_aff_add_pw_multi_aff(*upma, pma);
2971
2972         return 0;
2973 }
2974
2975 /* Return an isl_union_pw_multi_aff with the given "uset" as domain and
2976  * an unnamed zero-dimensional range.
2977  */
2978 __isl_give isl_union_pw_multi_aff *isl_union_pw_multi_aff_from_domain(
2979         __isl_take isl_union_set *uset)
2980 {
2981         isl_space *space;
2982         isl_union_pw_multi_aff *upma;
2983
2984         if (!uset)
2985                 return NULL;
2986
2987         space = isl_union_set_get_space(uset);
2988         upma = isl_union_pw_multi_aff_empty(space);
2989
2990         if (isl_union_set_foreach_set(uset,
2991                                     &add_pw_multi_aff_from_domain, &upma) < 0)
2992                 goto error;
2993
2994         isl_union_set_free(uset);
2995         return upma;
2996 error:
2997         isl_union_set_free(uset);
2998         isl_union_pw_multi_aff_free(upma);
2999         return NULL;
3000 }
3001
3002 /* Convert "pma" to an isl_map and add it to *umap.
3003  */
3004 static int map_from_pw_multi_aff(__isl_take isl_pw_multi_aff *pma, void *user)
3005 {
3006         isl_union_map **umap = user;
3007         isl_map *map;
3008
3009         map = isl_map_from_pw_multi_aff(pma);
3010         *umap = isl_union_map_add_map(*umap, map);
3011
3012         return 0;
3013 }
3014
3015 /* Construct a union map mapping the domain of the union
3016  * piecewise multi-affine expression to its range, with each dimension
3017  * in the range equated to the corresponding affine expression on its cell.
3018  */
3019 __isl_give isl_union_map *isl_union_map_from_union_pw_multi_aff(
3020         __isl_take isl_union_pw_multi_aff *upma)
3021 {
3022         isl_space *space;
3023         isl_union_map *umap;
3024
3025         if (!upma)
3026                 return NULL;
3027
3028         space = isl_union_pw_multi_aff_get_space(upma);
3029         umap = isl_union_map_empty(space);
3030
3031         if (isl_union_pw_multi_aff_foreach_pw_multi_aff(upma,
3032                                         &map_from_pw_multi_aff, &umap) < 0)
3033                 goto error;
3034
3035         isl_union_pw_multi_aff_free(upma);
3036         return umap;
3037 error:
3038         isl_union_pw_multi_aff_free(upma);
3039         isl_union_map_free(umap);
3040         return NULL;
3041 }
3042
3043 /* Local data for bin_entry and the callback "fn".
3044  */
3045 struct isl_union_pw_multi_aff_bin_data {
3046         isl_union_pw_multi_aff *upma2;
3047         isl_union_pw_multi_aff *res;
3048         isl_pw_multi_aff *pma;
3049         int (*fn)(void **entry, void *user);
3050 };
3051
3052 /* Given an isl_pw_multi_aff from upma1, store it in data->pma
3053  * and call data->fn for each isl_pw_multi_aff in data->upma2.
3054  */
3055 static int bin_entry(void **entry, void *user)
3056 {
3057         struct isl_union_pw_multi_aff_bin_data *data = user;
3058         isl_pw_multi_aff *pma = *entry;
3059
3060         data->pma = pma;
3061         if (isl_hash_table_foreach(data->upma2->dim->ctx, &data->upma2->table,
3062                                    data->fn, data) < 0)
3063                 return -1;
3064
3065         return 0;
3066 }
3067
3068 /* Call "fn" on each pair of isl_pw_multi_affs in "upma1" and "upma2".
3069  * The isl_pw_multi_aff from upma1 is stored in data->pma (where data is
3070  * passed as user field) and the isl_pw_multi_aff from upma2 is available
3071  * as *entry.  The callback should adjust data->res if desired.
3072  */
3073 static __isl_give isl_union_pw_multi_aff *bin_op(
3074         __isl_take isl_union_pw_multi_aff *upma1,
3075         __isl_take isl_union_pw_multi_aff *upma2,
3076         int (*fn)(void **entry, void *user))
3077 {
3078         isl_space *space;
3079         struct isl_union_pw_multi_aff_bin_data data = { NULL, NULL, NULL, fn };
3080
3081         space = isl_union_pw_multi_aff_get_space(upma2);
3082         upma1 = isl_union_pw_multi_aff_align_params(upma1, space);
3083         space = isl_union_pw_multi_aff_get_space(upma1);
3084         upma2 = isl_union_pw_multi_aff_align_params(upma2, space);
3085
3086         if (!upma1 || !upma2)
3087                 goto error;
3088
3089         data.upma2 = upma2;
3090         data.res = isl_union_pw_multi_aff_alloc(isl_space_copy(upma1->dim),
3091                                        upma1->table.n);
3092         if (isl_hash_table_foreach(upma1->dim->ctx, &upma1->table,
3093                                    &bin_entry, &data) < 0)
3094                 goto error;
3095
3096         isl_union_pw_multi_aff_free(upma1);
3097         isl_union_pw_multi_aff_free(upma2);
3098         return data.res;
3099 error:
3100         isl_union_pw_multi_aff_free(upma1);
3101         isl_union_pw_multi_aff_free(upma2);
3102         isl_union_pw_multi_aff_free(data.res);
3103         return NULL;
3104 }
3105
3106 /* Given two isl_multi_affs A -> B and C -> D,
3107  * construct an isl_multi_aff (A * C) -> (B, D).
3108  */
3109 __isl_give isl_multi_aff *isl_multi_aff_flat_range_product(
3110         __isl_take isl_multi_aff *ma1, __isl_take isl_multi_aff *ma2)
3111 {
3112         int i, n1, n2;
3113         isl_aff *aff;
3114         isl_space *space;
3115         isl_multi_aff *res;
3116
3117         if (!ma1 || !ma2)
3118                 goto error;
3119
3120         space = isl_space_range_product(isl_multi_aff_get_space(ma1),
3121                                         isl_multi_aff_get_space(ma2));
3122         space = isl_space_flatten_range(space);
3123         res = isl_multi_aff_alloc(space);
3124
3125         n1 = isl_multi_aff_dim(ma1, isl_dim_out);
3126         n2 = isl_multi_aff_dim(ma2, isl_dim_out);
3127
3128         for (i = 0; i < n1; ++i) {
3129                 aff = isl_multi_aff_get_aff(ma1, i);
3130                 res = isl_multi_aff_set_aff(res, i, aff);
3131         }
3132
3133         for (i = 0; i < n2; ++i) {
3134                 aff = isl_multi_aff_get_aff(ma2, i);
3135                 res = isl_multi_aff_set_aff(res, n1 + i, aff);
3136         }
3137
3138         isl_multi_aff_free(ma1);
3139         isl_multi_aff_free(ma2);
3140         return res;
3141 error:
3142         isl_multi_aff_free(ma1);
3143         isl_multi_aff_free(ma2);
3144         return NULL;
3145 }
3146
3147 /* Given two aligned isl_pw_multi_affs A -> B and C -> D,
3148  * construct an isl_pw_multi_aff (A * C) -> (B, D).
3149  */
3150 static __isl_give isl_pw_multi_aff *pw_multi_aff_flat_range_product(
3151         __isl_take isl_pw_multi_aff *pma1, __isl_take isl_pw_multi_aff *pma2)
3152 {
3153         isl_space *space;
3154
3155         space = isl_space_range_product(isl_pw_multi_aff_get_space(pma1),
3156                                         isl_pw_multi_aff_get_space(pma2));
3157         space = isl_space_flatten_range(space);
3158         return isl_pw_multi_aff_on_shared_domain_in(pma1, pma2, space,
3159                                             &isl_multi_aff_flat_range_product);
3160 }
3161
3162 /* Given two isl_pw_multi_affs A -> B and C -> D,
3163  * construct an isl_pw_multi_aff (A * C) -> (B, D).
3164  */
3165 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_flat_range_product(
3166         __isl_take isl_pw_multi_aff *pma1, __isl_take isl_pw_multi_aff *pma2)
3167 {
3168         return isl_pw_multi_aff_align_params_pw_pw_and(pma1, pma2,
3169                                             &pw_multi_aff_flat_range_product);
3170 }
3171
3172 /* If data->pma and *entry have the same domain space, then compute
3173  * their flat range product and the result to data->res.
3174  */
3175 static int flat_range_product_entry(void **entry, void *user)
3176 {
3177         struct isl_union_pw_multi_aff_bin_data *data = user;
3178         isl_pw_multi_aff *pma2 = *entry;
3179
3180         if (!isl_space_tuple_match(data->pma->dim, isl_dim_in,
3181                                  pma2->dim, isl_dim_in))
3182                 return 0;
3183
3184         pma2 = isl_pw_multi_aff_flat_range_product(
3185                                         isl_pw_multi_aff_copy(data->pma),
3186                                         isl_pw_multi_aff_copy(pma2));
3187
3188         data->res = isl_union_pw_multi_aff_add_pw_multi_aff(data->res, pma2);
3189
3190         return 0;
3191 }
3192
3193 /* Given two isl_union_pw_multi_affs A -> B and C -> D,
3194  * construct an isl_union_pw_multi_aff (A * C) -> (B, D).
3195  */
3196 __isl_give isl_union_pw_multi_aff *isl_union_pw_multi_aff_flat_range_product(
3197         __isl_take isl_union_pw_multi_aff *upma1,
3198         __isl_take isl_union_pw_multi_aff *upma2)
3199 {
3200         return bin_op(upma1, upma2, &flat_range_product_entry);
3201 }