isl_aff_normalize: sort divs in local space of affine expression
[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 MIT 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 /* Swap divs "a" and "b" in "aff", which is assumed to be non-NULL.
652  *
653  * Even though this function is only called on isl_affs with a single
654  * reference, we are careful to only change aff->v and aff->ls together.
655  */
656 static __isl_give isl_aff *swap_div(__isl_take isl_aff *aff, int a, int b)
657 {
658         unsigned off = isl_local_space_offset(aff->ls, isl_dim_div);
659         isl_local_space *ls;
660         isl_vec *v;
661
662         ls = isl_local_space_copy(aff->ls);
663         ls = isl_local_space_swap_div(ls, a, b);
664         v = isl_vec_copy(aff->v);
665         v = isl_vec_cow(v);
666         if (!ls || !v)
667                 goto error;
668
669         isl_int_swap(v->el[1 + off + a], v->el[1 + off + b]);
670         isl_vec_free(aff->v);
671         aff->v = v;
672         isl_local_space_free(aff->ls);
673         aff->ls = ls;
674
675         return aff;
676 error:
677         isl_vec_free(v);
678         isl_local_space_free(ls);
679         return isl_aff_free(aff);
680 }
681
682 /* Sort the divs in the local space of "aff" according to
683  * the comparison function "cmp_row" in isl_local_space.c
684  *
685  * Reordering divs does not change the semantics of "aff",
686  * so there is no need to call isl_aff_cow.
687  * Moreover, this function is currently only called on isl_affs
688  * with a single reference.
689  */
690 static __isl_give isl_aff *sort_divs(__isl_take isl_aff *aff)
691 {
692         int i, j, n;
693
694         if (!aff)
695                 return NULL;
696
697         n = isl_aff_dim(aff, isl_dim_div);
698         for (i = 1; i < n; ++i) {
699                 for (j = i - 1; j >= 0; --j) {
700                         if (isl_mat_cmp_div(aff->ls->div, j, j + 1) <= 0)
701                                 break;
702                         aff = swap_div(aff, j, j + 1);
703                         if (!aff)
704                                 return NULL;
705                 }
706         }
707
708         return aff;
709 }
710
711 /* Normalize the representation of "aff".
712  *
713  * This function should only be called of "new" isl_affs, i.e.,
714  * with only a single reference.  We therefore do not need to
715  * worry about affecting other instances.
716  */
717 __isl_give isl_aff *isl_aff_normalize(__isl_take isl_aff *aff)
718 {
719         if (!aff)
720                 return NULL;
721         aff->v = isl_vec_normalize(aff->v);
722         if (!aff->v)
723                 return isl_aff_free(aff);
724         aff = sort_divs(aff);
725         aff = isl_aff_remove_unused_divs(aff);
726         return aff;
727 }
728
729 /* Given f, return floor(f).
730  * If f is an integer expression, then just return f.
731  * If f is a constant, then return the constant floor(f).
732  * Otherwise, if f = g/m, write g = q m + r,
733  * create a new div d = [r/m] and return the expression q + d.
734  * The coefficients in r are taken to lie between -m/2 and m/2.
735  */
736 __isl_give isl_aff *isl_aff_floor(__isl_take isl_aff *aff)
737 {
738         int i;
739         int size;
740         isl_ctx *ctx;
741         isl_vec *div;
742
743         if (!aff)
744                 return NULL;
745
746         if (isl_int_is_one(aff->v->el[0]))
747                 return aff;
748
749         aff = isl_aff_cow(aff);
750         if (!aff)
751                 return NULL;
752
753         aff->v = isl_vec_cow(aff->v);
754         if (!aff->v)
755                 return isl_aff_free(aff);
756
757         if (isl_aff_is_cst(aff)) {
758                 isl_int_fdiv_q(aff->v->el[1], aff->v->el[1], aff->v->el[0]);
759                 isl_int_set_si(aff->v->el[0], 1);
760                 return aff;
761         }
762
763         div = isl_vec_copy(aff->v);
764         div = isl_vec_cow(div);
765         if (!div)
766                 return isl_aff_free(aff);
767
768         ctx = isl_aff_get_ctx(aff);
769         isl_int_fdiv_q(aff->v->el[0], aff->v->el[0], ctx->two);
770         for (i = 1; i < aff->v->size; ++i) {
771                 isl_int_fdiv_r(div->el[i], div->el[i], div->el[0]);
772                 isl_int_fdiv_q(aff->v->el[i], aff->v->el[i], div->el[0]);
773                 if (isl_int_gt(div->el[i], aff->v->el[0])) {
774                         isl_int_sub(div->el[i], div->el[i], div->el[0]);
775                         isl_int_add_ui(aff->v->el[i], aff->v->el[i], 1);
776                 }
777         }
778
779         aff->ls = isl_local_space_add_div(aff->ls, div);
780         if (!aff->ls)
781                 return isl_aff_free(aff);
782
783         size = aff->v->size;
784         aff->v = isl_vec_extend(aff->v, size + 1);
785         if (!aff->v)
786                 return isl_aff_free(aff);
787         isl_int_set_si(aff->v->el[0], 1);
788         isl_int_set_si(aff->v->el[size], 1);
789
790         return aff;
791 }
792
793 /* Compute
794  *
795  *      aff mod m = aff - m * floor(aff/m)
796  */
797 __isl_give isl_aff *isl_aff_mod(__isl_take isl_aff *aff, isl_int m)
798 {
799         isl_aff *res;
800
801         res = isl_aff_copy(aff);
802         aff = isl_aff_scale_down(aff, m);
803         aff = isl_aff_floor(aff);
804         aff = isl_aff_scale(aff, m);
805         res = isl_aff_sub(res, aff);
806
807         return res;
808 }
809
810 /* Compute
811  *
812  *      pwaff mod m = pwaff - m * floor(pwaff/m)
813  */
814 __isl_give isl_pw_aff *isl_pw_aff_mod(__isl_take isl_pw_aff *pwaff, isl_int m)
815 {
816         isl_pw_aff *res;
817
818         res = isl_pw_aff_copy(pwaff);
819         pwaff = isl_pw_aff_scale_down(pwaff, m);
820         pwaff = isl_pw_aff_floor(pwaff);
821         pwaff = isl_pw_aff_scale(pwaff, m);
822         res = isl_pw_aff_sub(res, pwaff);
823
824         return res;
825 }
826
827 /* Given f, return ceil(f).
828  * If f is an integer expression, then just return f.
829  * Otherwise, create a new div d = [-f] and return the expression -d.
830  */
831 __isl_give isl_aff *isl_aff_ceil(__isl_take isl_aff *aff)
832 {
833         if (!aff)
834                 return NULL;
835
836         if (isl_int_is_one(aff->v->el[0]))
837                 return aff;
838
839         aff = isl_aff_neg(aff);
840         aff = isl_aff_floor(aff);
841         aff = isl_aff_neg(aff);
842
843         return aff;
844 }
845
846 /* Apply the expansion computed by isl_merge_divs.
847  * The expansion itself is given by "exp" while the resulting
848  * list of divs is given by "div".
849  */
850 __isl_give isl_aff *isl_aff_expand_divs( __isl_take isl_aff *aff,
851         __isl_take isl_mat *div, int *exp)
852 {
853         int i, j;
854         int old_n_div;
855         int new_n_div;
856         int offset;
857
858         aff = isl_aff_cow(aff);
859         if (!aff || !div)
860                 goto error;
861
862         old_n_div = isl_local_space_dim(aff->ls, isl_dim_div);
863         new_n_div = isl_mat_rows(div);
864         if (new_n_div < old_n_div)
865                 isl_die(isl_mat_get_ctx(div), isl_error_invalid,
866                         "not an expansion", goto error);
867
868         aff->v = isl_vec_extend(aff->v, aff->v->size + new_n_div - old_n_div);
869         if (!aff->v)
870                 goto error;
871
872         offset = 1 + isl_local_space_offset(aff->ls, isl_dim_div);
873         j = old_n_div - 1;
874         for (i = new_n_div - 1; i >= 0; --i) {
875                 if (j >= 0 && exp[j] == i) {
876                         if (i != j)
877                                 isl_int_swap(aff->v->el[offset + i],
878                                              aff->v->el[offset + j]);
879                         j--;
880                 } else
881                         isl_int_set_si(aff->v->el[offset + i], 0);
882         }
883
884         aff->ls = isl_local_space_replace_divs(aff->ls, isl_mat_copy(div));
885         if (!aff->ls)
886                 goto error;
887         isl_mat_free(div);
888         return aff;
889 error:
890         isl_aff_free(aff);
891         isl_mat_free(div);
892         return NULL;
893 }
894
895 /* Add two affine expressions that live in the same local space.
896  */
897 static __isl_give isl_aff *add_expanded(__isl_take isl_aff *aff1,
898         __isl_take isl_aff *aff2)
899 {
900         isl_int gcd, f;
901
902         aff1 = isl_aff_cow(aff1);
903         if (!aff1 || !aff2)
904                 goto error;
905
906         aff1->v = isl_vec_cow(aff1->v);
907         if (!aff1->v)
908                 goto error;
909
910         isl_int_init(gcd);
911         isl_int_init(f);
912         isl_int_gcd(gcd, aff1->v->el[0], aff2->v->el[0]);
913         isl_int_divexact(f, aff2->v->el[0], gcd);
914         isl_seq_scale(aff1->v->el + 1, aff1->v->el + 1, f, aff1->v->size - 1);
915         isl_int_divexact(f, aff1->v->el[0], gcd);
916         isl_seq_addmul(aff1->v->el + 1, f, aff2->v->el + 1, aff1->v->size - 1);
917         isl_int_divexact(f, aff2->v->el[0], gcd);
918         isl_int_mul(aff1->v->el[0], aff1->v->el[0], f);
919         isl_int_clear(f);
920         isl_int_clear(gcd);
921
922         isl_aff_free(aff2);
923         return aff1;
924 error:
925         isl_aff_free(aff1);
926         isl_aff_free(aff2);
927         return NULL;
928 }
929
930 __isl_give isl_aff *isl_aff_add(__isl_take isl_aff *aff1,
931         __isl_take isl_aff *aff2)
932 {
933         isl_ctx *ctx;
934         int *exp1 = NULL;
935         int *exp2 = NULL;
936         isl_mat *div;
937
938         if (!aff1 || !aff2)
939                 goto error;
940
941         ctx = isl_aff_get_ctx(aff1);
942         if (!isl_space_is_equal(aff1->ls->dim, aff2->ls->dim))
943                 isl_die(ctx, isl_error_invalid,
944                         "spaces don't match", goto error);
945
946         if (aff1->ls->div->n_row == 0 && aff2->ls->div->n_row == 0)
947                 return add_expanded(aff1, aff2);
948
949         exp1 = isl_alloc_array(ctx, int, aff1->ls->div->n_row);
950         exp2 = isl_alloc_array(ctx, int, aff2->ls->div->n_row);
951         if (!exp1 || !exp2)
952                 goto error;
953
954         div = isl_merge_divs(aff1->ls->div, aff2->ls->div, exp1, exp2);
955         aff1 = isl_aff_expand_divs(aff1, isl_mat_copy(div), exp1);
956         aff2 = isl_aff_expand_divs(aff2, div, exp2);
957         free(exp1);
958         free(exp2);
959
960         return add_expanded(aff1, aff2);
961 error:
962         free(exp1);
963         free(exp2);
964         isl_aff_free(aff1);
965         isl_aff_free(aff2);
966         return NULL;
967 }
968
969 __isl_give isl_aff *isl_aff_sub(__isl_take isl_aff *aff1,
970         __isl_take isl_aff *aff2)
971 {
972         return isl_aff_add(aff1, isl_aff_neg(aff2));
973 }
974
975 __isl_give isl_aff *isl_aff_scale(__isl_take isl_aff *aff, isl_int f)
976 {
977         isl_int gcd;
978
979         if (isl_int_is_one(f))
980                 return aff;
981
982         aff = isl_aff_cow(aff);
983         if (!aff)
984                 return NULL;
985         aff->v = isl_vec_cow(aff->v);
986         if (!aff->v)
987                 return isl_aff_free(aff);
988
989         isl_int_init(gcd);
990         isl_int_gcd(gcd, aff->v->el[0], f);
991         isl_int_divexact(aff->v->el[0], aff->v->el[0], gcd);
992         isl_int_divexact(gcd, f, gcd);
993         isl_seq_scale(aff->v->el + 1, aff->v->el + 1, gcd, aff->v->size - 1);
994         isl_int_clear(gcd);
995
996         return aff;
997 }
998
999 __isl_give isl_aff *isl_aff_scale_down(__isl_take isl_aff *aff, isl_int f)
1000 {
1001         isl_int gcd;
1002
1003         if (isl_int_is_one(f))
1004                 return aff;
1005
1006         aff = isl_aff_cow(aff);
1007         if (!aff)
1008                 return NULL;
1009
1010         if (isl_int_is_zero(f))
1011                 isl_die(isl_aff_get_ctx(aff), isl_error_invalid,
1012                         "cannot scale down by zero", return isl_aff_free(aff));
1013
1014         aff->v = isl_vec_cow(aff->v);
1015         if (!aff->v)
1016                 return isl_aff_free(aff);
1017
1018         isl_int_init(gcd);
1019         isl_seq_gcd(aff->v->el + 1, aff->v->size - 1, &gcd);
1020         isl_int_gcd(gcd, gcd, f);
1021         isl_seq_scale_down(aff->v->el + 1, aff->v->el + 1, gcd, aff->v->size - 1);
1022         isl_int_divexact(gcd, f, gcd);
1023         isl_int_mul(aff->v->el[0], aff->v->el[0], gcd);
1024         isl_int_clear(gcd);
1025
1026         return aff;
1027 }
1028
1029 __isl_give isl_aff *isl_aff_scale_down_ui(__isl_take isl_aff *aff, unsigned f)
1030 {
1031         isl_int v;
1032
1033         if (f == 1)
1034                 return aff;
1035
1036         isl_int_init(v);
1037         isl_int_set_ui(v, f);
1038         aff = isl_aff_scale_down(aff, v);
1039         isl_int_clear(v);
1040
1041         return aff;
1042 }
1043
1044 __isl_give isl_aff *isl_aff_set_dim_name(__isl_take isl_aff *aff,
1045         enum isl_dim_type type, unsigned pos, const char *s)
1046 {
1047         aff = isl_aff_cow(aff);
1048         if (!aff)
1049                 return NULL;
1050         if (type == isl_dim_out)
1051                 isl_die(aff->v->ctx, isl_error_invalid,
1052                         "cannot set name of output/set dimension",
1053                         return isl_aff_free(aff));
1054         if (type == isl_dim_in)
1055                 type = isl_dim_set;
1056         aff->ls = isl_local_space_set_dim_name(aff->ls, type, pos, s);
1057         if (!aff->ls)
1058                 return isl_aff_free(aff);
1059
1060         return aff;
1061 }
1062
1063 __isl_give isl_aff *isl_aff_set_dim_id(__isl_take isl_aff *aff,
1064         enum isl_dim_type type, unsigned pos, __isl_take isl_id *id)
1065 {
1066         aff = isl_aff_cow(aff);
1067         if (!aff)
1068                 return isl_id_free(id);
1069         if (type == isl_dim_out)
1070                 isl_die(aff->v->ctx, isl_error_invalid,
1071                         "cannot set name of output/set dimension",
1072                         goto error);
1073         if (type == isl_dim_in)
1074                 type = isl_dim_set;
1075         aff->ls = isl_local_space_set_dim_id(aff->ls, type, pos, id);
1076         if (!aff->ls)
1077                 return isl_aff_free(aff);
1078
1079         return aff;
1080 error:
1081         isl_id_free(id);
1082         isl_aff_free(aff);
1083         return NULL;
1084 }
1085
1086 /* Exploit the equalities in "eq" to simplify the affine expression
1087  * and the expressions of the integer divisions in the local space.
1088  * The integer divisions in this local space are assumed to appear
1089  * as regular dimensions in "eq".
1090  */
1091 static __isl_give isl_aff *isl_aff_substitute_equalities_lifted(
1092         __isl_take isl_aff *aff, __isl_take isl_basic_set *eq)
1093 {
1094         int i, j;
1095         unsigned total;
1096         unsigned n_div;
1097
1098         if (!eq)
1099                 goto error;
1100         if (eq->n_eq == 0) {
1101                 isl_basic_set_free(eq);
1102                 return aff;
1103         }
1104
1105         aff = isl_aff_cow(aff);
1106         if (!aff)
1107                 goto error;
1108
1109         aff->ls = isl_local_space_substitute_equalities(aff->ls,
1110                                                         isl_basic_set_copy(eq));
1111         if (!aff->ls)
1112                 goto error;
1113
1114         total = 1 + isl_space_dim(eq->dim, isl_dim_all);
1115         n_div = eq->n_div;
1116         for (i = 0; i < eq->n_eq; ++i) {
1117                 j = isl_seq_last_non_zero(eq->eq[i], total + n_div);
1118                 if (j < 0 || j == 0 || j >= total)
1119                         continue;
1120
1121                 isl_seq_elim(aff->v->el + 1, eq->eq[i], j, total,
1122                                 &aff->v->el[0]);
1123         }
1124
1125         isl_basic_set_free(eq);
1126         aff = isl_aff_normalize(aff);
1127         return aff;
1128 error:
1129         isl_basic_set_free(eq);
1130         isl_aff_free(aff);
1131         return NULL;
1132 }
1133
1134 /* Exploit the equalities in "eq" to simplify the affine expression
1135  * and the expressions of the integer divisions in the local space.
1136  */
1137 static __isl_give isl_aff *isl_aff_substitute_equalities(
1138         __isl_take isl_aff *aff, __isl_take isl_basic_set *eq)
1139 {
1140         int n_div;
1141
1142         if (!aff || !eq)
1143                 goto error;
1144         n_div = isl_local_space_dim(aff->ls, isl_dim_div);
1145         if (n_div > 0)
1146                 eq = isl_basic_set_add(eq, isl_dim_set, n_div);
1147         return isl_aff_substitute_equalities_lifted(aff, eq);
1148 error:
1149         isl_basic_set_free(eq);
1150         isl_aff_free(aff);
1151         return NULL;
1152 }
1153
1154 /* Look for equalities among the variables shared by context and aff
1155  * and the integer divisions of aff, if any.
1156  * The equalities are then used to eliminate coefficients and/or integer
1157  * divisions from aff.
1158  */
1159 __isl_give isl_aff *isl_aff_gist(__isl_take isl_aff *aff,
1160         __isl_take isl_set *context)
1161 {
1162         isl_basic_set *hull;
1163         int n_div;
1164
1165         if (!aff)
1166                 goto error;
1167         n_div = isl_local_space_dim(aff->ls, isl_dim_div);
1168         if (n_div > 0) {
1169                 isl_basic_set *bset;
1170                 isl_local_space *ls;
1171                 context = isl_set_add_dims(context, isl_dim_set, n_div);
1172                 ls = isl_aff_get_domain_local_space(aff);
1173                 bset = isl_basic_set_from_local_space(ls);
1174                 bset = isl_basic_set_lift(bset);
1175                 bset = isl_basic_set_flatten(bset);
1176                 context = isl_set_intersect(context,
1177                                             isl_set_from_basic_set(bset));
1178         }
1179
1180         hull = isl_set_affine_hull(context);
1181         return isl_aff_substitute_equalities_lifted(aff, hull);
1182 error:
1183         isl_aff_free(aff);
1184         isl_set_free(context);
1185         return NULL;
1186 }
1187
1188 __isl_give isl_aff *isl_aff_gist_params(__isl_take isl_aff *aff,
1189         __isl_take isl_set *context)
1190 {
1191         isl_set *dom_context = isl_set_universe(isl_aff_get_domain_space(aff));
1192         dom_context = isl_set_intersect_params(dom_context, context);
1193         return isl_aff_gist(aff, dom_context);
1194 }
1195
1196 /* Return a basic set containing those elements in the space
1197  * of aff where it is non-negative.
1198  */
1199 __isl_give isl_basic_set *isl_aff_nonneg_basic_set(__isl_take isl_aff *aff)
1200 {
1201         isl_constraint *ineq;
1202         isl_basic_set *bset;
1203
1204         ineq = isl_inequality_from_aff(aff);
1205
1206         bset = isl_basic_set_from_constraint(ineq);
1207         bset = isl_basic_set_simplify(bset);
1208         return bset;
1209 }
1210
1211 /* Return a basic set containing those elements in the domain space
1212  * of aff where it is negative.
1213  */
1214 __isl_give isl_basic_set *isl_aff_neg_basic_set(__isl_take isl_aff *aff)
1215 {
1216         aff = isl_aff_neg(aff);
1217         aff = isl_aff_add_constant_num_si(aff, -1);
1218         return isl_aff_nonneg_basic_set(aff);
1219 }
1220
1221 /* Return a basic set containing those elements in the space
1222  * of aff where it is zero.
1223  */
1224 __isl_give isl_basic_set *isl_aff_zero_basic_set(__isl_take isl_aff *aff)
1225 {
1226         isl_constraint *ineq;
1227         isl_basic_set *bset;
1228
1229         ineq = isl_equality_from_aff(aff);
1230
1231         bset = isl_basic_set_from_constraint(ineq);
1232         bset = isl_basic_set_simplify(bset);
1233         return bset;
1234 }
1235
1236 /* Return a basic set containing those elements in the shared space
1237  * of aff1 and aff2 where aff1 is greater than or equal to aff2.
1238  */
1239 __isl_give isl_basic_set *isl_aff_ge_basic_set(__isl_take isl_aff *aff1,
1240         __isl_take isl_aff *aff2)
1241 {
1242         aff1 = isl_aff_sub(aff1, aff2);
1243
1244         return isl_aff_nonneg_basic_set(aff1);
1245 }
1246
1247 /* Return a basic set containing those elements in the shared space
1248  * of aff1 and aff2 where aff1 is smaller than or equal to aff2.
1249  */
1250 __isl_give isl_basic_set *isl_aff_le_basic_set(__isl_take isl_aff *aff1,
1251         __isl_take isl_aff *aff2)
1252 {
1253         return isl_aff_ge_basic_set(aff2, aff1);
1254 }
1255
1256 __isl_give isl_aff *isl_aff_add_on_domain(__isl_keep isl_set *dom,
1257         __isl_take isl_aff *aff1, __isl_take isl_aff *aff2)
1258 {
1259         aff1 = isl_aff_add(aff1, aff2);
1260         aff1 = isl_aff_gist(aff1, isl_set_copy(dom));
1261         return aff1;
1262 }
1263
1264 int isl_aff_is_empty(__isl_keep isl_aff *aff)
1265 {
1266         if (!aff)
1267                 return -1;
1268
1269         return 0;
1270 }
1271
1272 /* Check whether the given affine expression has non-zero coefficient
1273  * for any dimension in the given range or if any of these dimensions
1274  * appear with non-zero coefficients in any of the integer divisions
1275  * involved in the affine expression.
1276  */
1277 int isl_aff_involves_dims(__isl_keep isl_aff *aff,
1278         enum isl_dim_type type, unsigned first, unsigned n)
1279 {
1280         int i;
1281         isl_ctx *ctx;
1282         int *active = NULL;
1283         int involves = 0;
1284
1285         if (!aff)
1286                 return -1;
1287         if (n == 0)
1288                 return 0;
1289
1290         ctx = isl_aff_get_ctx(aff);
1291         if (first + n > isl_aff_dim(aff, type))
1292                 isl_die(ctx, isl_error_invalid,
1293                         "range out of bounds", return -1);
1294
1295         active = isl_local_space_get_active(aff->ls, aff->v->el + 2);
1296         if (!active)
1297                 goto error;
1298
1299         first += isl_local_space_offset(aff->ls, type) - 1;
1300         for (i = 0; i < n; ++i)
1301                 if (active[first + i]) {
1302                         involves = 1;
1303                         break;
1304                 }
1305
1306         free(active);
1307
1308         return involves;
1309 error:
1310         free(active);
1311         return -1;
1312 }
1313
1314 __isl_give isl_aff *isl_aff_drop_dims(__isl_take isl_aff *aff,
1315         enum isl_dim_type type, unsigned first, unsigned n)
1316 {
1317         isl_ctx *ctx;
1318
1319         if (!aff)
1320                 return NULL;
1321         if (type == isl_dim_out)
1322                 isl_die(aff->v->ctx, isl_error_invalid,
1323                         "cannot drop output/set dimension",
1324                         return isl_aff_free(aff));
1325         if (type == isl_dim_in)
1326                 type = isl_dim_set;
1327         if (n == 0 && !isl_local_space_is_named_or_nested(aff->ls, type))
1328                 return aff;
1329
1330         ctx = isl_aff_get_ctx(aff);
1331         if (first + n > isl_local_space_dim(aff->ls, type))
1332                 isl_die(ctx, isl_error_invalid, "range out of bounds",
1333                         return isl_aff_free(aff));
1334
1335         aff = isl_aff_cow(aff);
1336         if (!aff)
1337                 return NULL;
1338
1339         aff->ls = isl_local_space_drop_dims(aff->ls, type, first, n);
1340         if (!aff->ls)
1341                 return isl_aff_free(aff);
1342
1343         first += 1 + isl_local_space_offset(aff->ls, type);
1344         aff->v = isl_vec_drop_els(aff->v, first, n);
1345         if (!aff->v)
1346                 return isl_aff_free(aff);
1347
1348         return aff;
1349 }
1350
1351 /* Project the domain of the affine expression onto its parameter space.
1352  * The affine expression may not involve any of the domain dimensions.
1353  */
1354 __isl_give isl_aff *isl_aff_project_domain_on_params(__isl_take isl_aff *aff)
1355 {
1356         isl_space *space;
1357         unsigned n;
1358         int involves;
1359
1360         n = isl_aff_dim(aff, isl_dim_in);
1361         involves = isl_aff_involves_dims(aff, isl_dim_in, 0, n);
1362         if (involves < 0)
1363                 return isl_aff_free(aff);
1364         if (involves)
1365                 isl_die(isl_aff_get_ctx(aff), isl_error_invalid,
1366                     "affine expression involves some of the domain dimensions",
1367                     return isl_aff_free(aff));
1368         aff = isl_aff_drop_dims(aff, isl_dim_in, 0, n);
1369         space = isl_aff_get_domain_space(aff);
1370         space = isl_space_params(space);
1371         aff = isl_aff_reset_domain_space(aff, space);
1372         return aff;
1373 }
1374
1375 __isl_give isl_aff *isl_aff_insert_dims(__isl_take isl_aff *aff,
1376         enum isl_dim_type type, unsigned first, unsigned n)
1377 {
1378         isl_ctx *ctx;
1379
1380         if (!aff)
1381                 return NULL;
1382         if (type == isl_dim_out)
1383                 isl_die(aff->v->ctx, isl_error_invalid,
1384                         "cannot insert output/set dimensions",
1385                         return isl_aff_free(aff));
1386         if (type == isl_dim_in)
1387                 type = isl_dim_set;
1388         if (n == 0 && !isl_local_space_is_named_or_nested(aff->ls, type))
1389                 return aff;
1390
1391         ctx = isl_aff_get_ctx(aff);
1392         if (first > isl_local_space_dim(aff->ls, type))
1393                 isl_die(ctx, isl_error_invalid, "position out of bounds",
1394                         return isl_aff_free(aff));
1395
1396         aff = isl_aff_cow(aff);
1397         if (!aff)
1398                 return NULL;
1399
1400         aff->ls = isl_local_space_insert_dims(aff->ls, type, first, n);
1401         if (!aff->ls)
1402                 return isl_aff_free(aff);
1403
1404         first += 1 + isl_local_space_offset(aff->ls, type);
1405         aff->v = isl_vec_insert_zero_els(aff->v, first, n);
1406         if (!aff->v)
1407                 return isl_aff_free(aff);
1408
1409         return aff;
1410 }
1411
1412 __isl_give isl_aff *isl_aff_add_dims(__isl_take isl_aff *aff,
1413         enum isl_dim_type type, unsigned n)
1414 {
1415         unsigned pos;
1416
1417         pos = isl_aff_dim(aff, type);
1418
1419         return isl_aff_insert_dims(aff, type, pos, n);
1420 }
1421
1422 __isl_give isl_pw_aff *isl_pw_aff_add_dims(__isl_take isl_pw_aff *pwaff,
1423         enum isl_dim_type type, unsigned n)
1424 {
1425         unsigned pos;
1426
1427         pos = isl_pw_aff_dim(pwaff, type);
1428
1429         return isl_pw_aff_insert_dims(pwaff, type, pos, n);
1430 }
1431
1432 __isl_give isl_pw_aff *isl_pw_aff_from_aff(__isl_take isl_aff *aff)
1433 {
1434         isl_set *dom = isl_set_universe(isl_aff_get_domain_space(aff));
1435         return isl_pw_aff_alloc(dom, aff);
1436 }
1437
1438 #undef PW
1439 #define PW isl_pw_aff
1440 #undef EL
1441 #define EL isl_aff
1442 #undef EL_IS_ZERO
1443 #define EL_IS_ZERO is_empty
1444 #undef ZERO
1445 #define ZERO empty
1446 #undef IS_ZERO
1447 #define IS_ZERO is_empty
1448 #undef FIELD
1449 #define FIELD aff
1450 #undef DEFAULT_IS_ZERO
1451 #define DEFAULT_IS_ZERO 0
1452
1453 #define NO_EVAL
1454 #define NO_OPT
1455 #define NO_MOVE_DIMS
1456 #define NO_LIFT
1457 #define NO_MORPH
1458
1459 #include <isl_pw_templ.c>
1460
1461 static __isl_give isl_set *align_params_pw_pw_set_and(
1462         __isl_take isl_pw_aff *pwaff1, __isl_take isl_pw_aff *pwaff2,
1463         __isl_give isl_set *(*fn)(__isl_take isl_pw_aff *pwaff1,
1464                                     __isl_take isl_pw_aff *pwaff2))
1465 {
1466         if (!pwaff1 || !pwaff2)
1467                 goto error;
1468         if (isl_space_match(pwaff1->dim, isl_dim_param,
1469                           pwaff2->dim, isl_dim_param))
1470                 return fn(pwaff1, pwaff2);
1471         if (!isl_space_has_named_params(pwaff1->dim) ||
1472             !isl_space_has_named_params(pwaff2->dim))
1473                 isl_die(isl_pw_aff_get_ctx(pwaff1), isl_error_invalid,
1474                         "unaligned unnamed parameters", goto error);
1475         pwaff1 = isl_pw_aff_align_params(pwaff1, isl_pw_aff_get_space(pwaff2));
1476         pwaff2 = isl_pw_aff_align_params(pwaff2, isl_pw_aff_get_space(pwaff1));
1477         return fn(pwaff1, pwaff2);
1478 error:
1479         isl_pw_aff_free(pwaff1);
1480         isl_pw_aff_free(pwaff2);
1481         return NULL;
1482 }
1483
1484 /* Compute a piecewise quasi-affine expression with a domain that
1485  * is the union of those of pwaff1 and pwaff2 and such that on each
1486  * cell, the quasi-affine expression is the better (according to cmp)
1487  * of those of pwaff1 and pwaff2.  If only one of pwaff1 or pwaff2
1488  * is defined on a given cell, then the associated expression
1489  * is the defined one.
1490  */
1491 static __isl_give isl_pw_aff *pw_aff_union_opt(__isl_take isl_pw_aff *pwaff1,
1492         __isl_take isl_pw_aff *pwaff2,
1493         __isl_give isl_basic_set *(*cmp)(__isl_take isl_aff *aff1,
1494                                         __isl_take isl_aff *aff2))
1495 {
1496         int i, j, n;
1497         isl_pw_aff *res;
1498         isl_ctx *ctx;
1499         isl_set *set;
1500
1501         if (!pwaff1 || !pwaff2)
1502                 goto error;
1503
1504         ctx = isl_space_get_ctx(pwaff1->dim);
1505         if (!isl_space_is_equal(pwaff1->dim, pwaff2->dim))
1506                 isl_die(ctx, isl_error_invalid,
1507                         "arguments should live in same space", goto error);
1508
1509         if (isl_pw_aff_is_empty(pwaff1)) {
1510                 isl_pw_aff_free(pwaff1);
1511                 return pwaff2;
1512         }
1513
1514         if (isl_pw_aff_is_empty(pwaff2)) {
1515                 isl_pw_aff_free(pwaff2);
1516                 return pwaff1;
1517         }
1518
1519         n = 2 * (pwaff1->n + 1) * (pwaff2->n + 1);
1520         res = isl_pw_aff_alloc_size(isl_space_copy(pwaff1->dim), n);
1521
1522         for (i = 0; i < pwaff1->n; ++i) {
1523                 set = isl_set_copy(pwaff1->p[i].set);
1524                 for (j = 0; j < pwaff2->n; ++j) {
1525                         struct isl_set *common;
1526                         isl_set *better;
1527
1528                         common = isl_set_intersect(
1529                                         isl_set_copy(pwaff1->p[i].set),
1530                                         isl_set_copy(pwaff2->p[j].set));
1531                         better = isl_set_from_basic_set(cmp(
1532                                         isl_aff_copy(pwaff2->p[j].aff),
1533                                         isl_aff_copy(pwaff1->p[i].aff)));
1534                         better = isl_set_intersect(common, better);
1535                         if (isl_set_plain_is_empty(better)) {
1536                                 isl_set_free(better);
1537                                 continue;
1538                         }
1539                         set = isl_set_subtract(set, isl_set_copy(better));
1540
1541                         res = isl_pw_aff_add_piece(res, better,
1542                                                 isl_aff_copy(pwaff2->p[j].aff));
1543                 }
1544                 res = isl_pw_aff_add_piece(res, set,
1545                                                 isl_aff_copy(pwaff1->p[i].aff));
1546         }
1547
1548         for (j = 0; j < pwaff2->n; ++j) {
1549                 set = isl_set_copy(pwaff2->p[j].set);
1550                 for (i = 0; i < pwaff1->n; ++i)
1551                         set = isl_set_subtract(set,
1552                                         isl_set_copy(pwaff1->p[i].set));
1553                 res = isl_pw_aff_add_piece(res, set,
1554                                                 isl_aff_copy(pwaff2->p[j].aff));
1555         }
1556
1557         isl_pw_aff_free(pwaff1);
1558         isl_pw_aff_free(pwaff2);
1559
1560         return res;
1561 error:
1562         isl_pw_aff_free(pwaff1);
1563         isl_pw_aff_free(pwaff2);
1564         return NULL;
1565 }
1566
1567 /* Compute a piecewise quasi-affine expression with a domain that
1568  * is the union of those of pwaff1 and pwaff2 and such that on each
1569  * cell, the quasi-affine expression is the maximum of those of pwaff1
1570  * and pwaff2.  If only one of pwaff1 or pwaff2 is defined on a given
1571  * cell, then the associated expression is the defined one.
1572  */
1573 static __isl_give isl_pw_aff *pw_aff_union_max(__isl_take isl_pw_aff *pwaff1,
1574         __isl_take isl_pw_aff *pwaff2)
1575 {
1576         return pw_aff_union_opt(pwaff1, pwaff2, &isl_aff_ge_basic_set);
1577 }
1578
1579 __isl_give isl_pw_aff *isl_pw_aff_union_max(__isl_take isl_pw_aff *pwaff1,
1580         __isl_take isl_pw_aff *pwaff2)
1581 {
1582         return isl_pw_aff_align_params_pw_pw_and(pwaff1, pwaff2,
1583                                                         &pw_aff_union_max);
1584 }
1585
1586 /* Compute a piecewise quasi-affine expression with a domain that
1587  * is the union of those of pwaff1 and pwaff2 and such that on each
1588  * cell, the quasi-affine expression is the minimum of those of pwaff1
1589  * and pwaff2.  If only one of pwaff1 or pwaff2 is defined on a given
1590  * cell, then the associated expression is the defined one.
1591  */
1592 static __isl_give isl_pw_aff *pw_aff_union_min(__isl_take isl_pw_aff *pwaff1,
1593         __isl_take isl_pw_aff *pwaff2)
1594 {
1595         return pw_aff_union_opt(pwaff1, pwaff2, &isl_aff_le_basic_set);
1596 }
1597
1598 __isl_give isl_pw_aff *isl_pw_aff_union_min(__isl_take isl_pw_aff *pwaff1,
1599         __isl_take isl_pw_aff *pwaff2)
1600 {
1601         return isl_pw_aff_align_params_pw_pw_and(pwaff1, pwaff2,
1602                                                         &pw_aff_union_min);
1603 }
1604
1605 __isl_give isl_pw_aff *isl_pw_aff_union_opt(__isl_take isl_pw_aff *pwaff1,
1606         __isl_take isl_pw_aff *pwaff2, int max)
1607 {
1608         if (max)
1609                 return isl_pw_aff_union_max(pwaff1, pwaff2);
1610         else
1611                 return isl_pw_aff_union_min(pwaff1, pwaff2);
1612 }
1613
1614 /* Construct a map with as domain the domain of pwaff and
1615  * one-dimensional range corresponding to the affine expressions.
1616  */
1617 static __isl_give isl_map *map_from_pw_aff(__isl_take isl_pw_aff *pwaff)
1618 {
1619         int i;
1620         isl_space *dim;
1621         isl_map *map;
1622
1623         if (!pwaff)
1624                 return NULL;
1625
1626         dim = isl_pw_aff_get_space(pwaff);
1627         map = isl_map_empty(dim);
1628
1629         for (i = 0; i < pwaff->n; ++i) {
1630                 isl_basic_map *bmap;
1631                 isl_map *map_i;
1632
1633                 bmap = isl_basic_map_from_aff(isl_aff_copy(pwaff->p[i].aff));
1634                 map_i = isl_map_from_basic_map(bmap);
1635                 map_i = isl_map_intersect_domain(map_i,
1636                                                 isl_set_copy(pwaff->p[i].set));
1637                 map = isl_map_union_disjoint(map, map_i);
1638         }
1639
1640         isl_pw_aff_free(pwaff);
1641
1642         return map;
1643 }
1644
1645 /* Construct a map with as domain the domain of pwaff and
1646  * one-dimensional range corresponding to the affine expressions.
1647  */
1648 __isl_give isl_map *isl_map_from_pw_aff(__isl_take isl_pw_aff *pwaff)
1649 {
1650         if (!pwaff)
1651                 return NULL;
1652         if (isl_space_is_set(pwaff->dim))
1653                 isl_die(isl_pw_aff_get_ctx(pwaff), isl_error_invalid,
1654                         "space of input is not a map",
1655                         return isl_pw_aff_free(pwaff));
1656         return map_from_pw_aff(pwaff);
1657 }
1658
1659 /* Construct a one-dimensional set with as parameter domain
1660  * the domain of pwaff and the single set dimension
1661  * corresponding to the affine expressions.
1662  */
1663 __isl_give isl_set *isl_set_from_pw_aff(__isl_take isl_pw_aff *pwaff)
1664 {
1665         if (!pwaff)
1666                 return NULL;
1667         if (!isl_space_is_set(pwaff->dim))
1668                 isl_die(isl_pw_aff_get_ctx(pwaff), isl_error_invalid,
1669                         "space of input is not a set",
1670                         return isl_pw_aff_free(pwaff));
1671         return map_from_pw_aff(pwaff);
1672 }
1673
1674 /* Return a set containing those elements in the domain
1675  * of pwaff where it is non-negative.
1676  */
1677 __isl_give isl_set *isl_pw_aff_nonneg_set(__isl_take isl_pw_aff *pwaff)
1678 {
1679         int i;
1680         isl_set *set;
1681
1682         if (!pwaff)
1683                 return NULL;
1684
1685         set = isl_set_empty(isl_pw_aff_get_domain_space(pwaff));
1686
1687         for (i = 0; i < pwaff->n; ++i) {
1688                 isl_basic_set *bset;
1689                 isl_set *set_i;
1690
1691                 bset = isl_aff_nonneg_basic_set(isl_aff_copy(pwaff->p[i].aff));
1692                 set_i = isl_set_from_basic_set(bset);
1693                 set_i = isl_set_intersect(set_i, isl_set_copy(pwaff->p[i].set));
1694                 set = isl_set_union_disjoint(set, set_i);
1695         }
1696
1697         isl_pw_aff_free(pwaff);
1698
1699         return set;
1700 }
1701
1702 /* Return a set containing those elements in the domain
1703  * of pwaff where it is zero (if complement is 0) or not zero
1704  * (if complement is 1).
1705  */
1706 static __isl_give isl_set *pw_aff_zero_set(__isl_take isl_pw_aff *pwaff,
1707         int complement)
1708 {
1709         int i;
1710         isl_set *set;
1711
1712         if (!pwaff)
1713                 return NULL;
1714
1715         set = isl_set_empty(isl_pw_aff_get_domain_space(pwaff));
1716
1717         for (i = 0; i < pwaff->n; ++i) {
1718                 isl_basic_set *bset;
1719                 isl_set *set_i, *zero;
1720
1721                 bset = isl_aff_zero_basic_set(isl_aff_copy(pwaff->p[i].aff));
1722                 zero = isl_set_from_basic_set(bset);
1723                 set_i = isl_set_copy(pwaff->p[i].set);
1724                 if (complement)
1725                         set_i = isl_set_subtract(set_i, zero);
1726                 else
1727                         set_i = isl_set_intersect(set_i, zero);
1728                 set = isl_set_union_disjoint(set, set_i);
1729         }
1730
1731         isl_pw_aff_free(pwaff);
1732
1733         return set;
1734 }
1735
1736 /* Return a set containing those elements in the domain
1737  * of pwaff where it is zero.
1738  */
1739 __isl_give isl_set *isl_pw_aff_zero_set(__isl_take isl_pw_aff *pwaff)
1740 {
1741         return pw_aff_zero_set(pwaff, 0);
1742 }
1743
1744 /* Return a set containing those elements in the domain
1745  * of pwaff where it is not zero.
1746  */
1747 __isl_give isl_set *isl_pw_aff_non_zero_set(__isl_take isl_pw_aff *pwaff)
1748 {
1749         return pw_aff_zero_set(pwaff, 1);
1750 }
1751
1752 /* Return a set containing those elements in the shared domain
1753  * of pwaff1 and pwaff2 where pwaff1 is greater than (or equal) to pwaff2.
1754  *
1755  * We compute the difference on the shared domain and then construct
1756  * the set of values where this difference is non-negative.
1757  * If strict is set, we first subtract 1 from the difference.
1758  * If equal is set, we only return the elements where pwaff1 and pwaff2
1759  * are equal.
1760  */
1761 static __isl_give isl_set *pw_aff_gte_set(__isl_take isl_pw_aff *pwaff1,
1762         __isl_take isl_pw_aff *pwaff2, int strict, int equal)
1763 {
1764         isl_set *set1, *set2;
1765
1766         set1 = isl_pw_aff_domain(isl_pw_aff_copy(pwaff1));
1767         set2 = isl_pw_aff_domain(isl_pw_aff_copy(pwaff2));
1768         set1 = isl_set_intersect(set1, set2);
1769         pwaff1 = isl_pw_aff_intersect_domain(pwaff1, isl_set_copy(set1));
1770         pwaff2 = isl_pw_aff_intersect_domain(pwaff2, isl_set_copy(set1));
1771         pwaff1 = isl_pw_aff_add(pwaff1, isl_pw_aff_neg(pwaff2));
1772
1773         if (strict) {
1774                 isl_space *dim = isl_set_get_space(set1);
1775                 isl_aff *aff;
1776                 aff = isl_aff_zero_on_domain(isl_local_space_from_space(dim));
1777                 aff = isl_aff_add_constant_si(aff, -1);
1778                 pwaff1 = isl_pw_aff_add(pwaff1, isl_pw_aff_alloc(set1, aff));
1779         } else
1780                 isl_set_free(set1);
1781
1782         if (equal)
1783                 return isl_pw_aff_zero_set(pwaff1);
1784         return isl_pw_aff_nonneg_set(pwaff1);
1785 }
1786
1787 /* Return a set containing those elements in the shared domain
1788  * of pwaff1 and pwaff2 where pwaff1 is equal to pwaff2.
1789  */
1790 static __isl_give isl_set *pw_aff_eq_set(__isl_take isl_pw_aff *pwaff1,
1791         __isl_take isl_pw_aff *pwaff2)
1792 {
1793         return pw_aff_gte_set(pwaff1, pwaff2, 0, 1);
1794 }
1795
1796 __isl_give isl_set *isl_pw_aff_eq_set(__isl_take isl_pw_aff *pwaff1,
1797         __isl_take isl_pw_aff *pwaff2)
1798 {
1799         return align_params_pw_pw_set_and(pwaff1, pwaff2, &pw_aff_eq_set);
1800 }
1801
1802 /* Return a set containing those elements in the shared domain
1803  * of pwaff1 and pwaff2 where pwaff1 is greater than or equal to pwaff2.
1804  */
1805 static __isl_give isl_set *pw_aff_ge_set(__isl_take isl_pw_aff *pwaff1,
1806         __isl_take isl_pw_aff *pwaff2)
1807 {
1808         return pw_aff_gte_set(pwaff1, pwaff2, 0, 0);
1809 }
1810
1811 __isl_give isl_set *isl_pw_aff_ge_set(__isl_take isl_pw_aff *pwaff1,
1812         __isl_take isl_pw_aff *pwaff2)
1813 {
1814         return align_params_pw_pw_set_and(pwaff1, pwaff2, &pw_aff_ge_set);
1815 }
1816
1817 /* Return a set containing those elements in the shared domain
1818  * of pwaff1 and pwaff2 where pwaff1 is strictly greater than pwaff2.
1819  */
1820 static __isl_give isl_set *pw_aff_gt_set(__isl_take isl_pw_aff *pwaff1,
1821         __isl_take isl_pw_aff *pwaff2)
1822 {
1823         return pw_aff_gte_set(pwaff1, pwaff2, 1, 0);
1824 }
1825
1826 __isl_give isl_set *isl_pw_aff_gt_set(__isl_take isl_pw_aff *pwaff1,
1827         __isl_take isl_pw_aff *pwaff2)
1828 {
1829         return align_params_pw_pw_set_and(pwaff1, pwaff2, &pw_aff_gt_set);
1830 }
1831
1832 __isl_give isl_set *isl_pw_aff_le_set(__isl_take isl_pw_aff *pwaff1,
1833         __isl_take isl_pw_aff *pwaff2)
1834 {
1835         return isl_pw_aff_ge_set(pwaff2, pwaff1);
1836 }
1837
1838 __isl_give isl_set *isl_pw_aff_lt_set(__isl_take isl_pw_aff *pwaff1,
1839         __isl_take isl_pw_aff *pwaff2)
1840 {
1841         return isl_pw_aff_gt_set(pwaff2, pwaff1);
1842 }
1843
1844 /* Return a set containing those elements in the shared domain
1845  * of the elements of list1 and list2 where each element in list1
1846  * has the relation specified by "fn" with each element in list2.
1847  */
1848 static __isl_give isl_set *pw_aff_list_set(__isl_take isl_pw_aff_list *list1,
1849         __isl_take isl_pw_aff_list *list2,
1850         __isl_give isl_set *(*fn)(__isl_take isl_pw_aff *pwaff1,
1851                                     __isl_take isl_pw_aff *pwaff2))
1852 {
1853         int i, j;
1854         isl_ctx *ctx;
1855         isl_set *set;
1856
1857         if (!list1 || !list2)
1858                 goto error;
1859
1860         ctx = isl_pw_aff_list_get_ctx(list1);
1861         if (list1->n < 1 || list2->n < 1)
1862                 isl_die(ctx, isl_error_invalid,
1863                         "list should contain at least one element", goto error);
1864
1865         set = isl_set_universe(isl_pw_aff_get_domain_space(list1->p[0]));
1866         for (i = 0; i < list1->n; ++i)
1867                 for (j = 0; j < list2->n; ++j) {
1868                         isl_set *set_ij;
1869
1870                         set_ij = fn(isl_pw_aff_copy(list1->p[i]),
1871                                     isl_pw_aff_copy(list2->p[j]));
1872                         set = isl_set_intersect(set, set_ij);
1873                 }
1874
1875         isl_pw_aff_list_free(list1);
1876         isl_pw_aff_list_free(list2);
1877         return set;
1878 error:
1879         isl_pw_aff_list_free(list1);
1880         isl_pw_aff_list_free(list2);
1881         return NULL;
1882 }
1883
1884 /* Return a set containing those elements in the shared domain
1885  * of the elements of list1 and list2 where each element in list1
1886  * is equal to each element in list2.
1887  */
1888 __isl_give isl_set *isl_pw_aff_list_eq_set(__isl_take isl_pw_aff_list *list1,
1889         __isl_take isl_pw_aff_list *list2)
1890 {
1891         return pw_aff_list_set(list1, list2, &isl_pw_aff_eq_set);
1892 }
1893
1894 __isl_give isl_set *isl_pw_aff_list_ne_set(__isl_take isl_pw_aff_list *list1,
1895         __isl_take isl_pw_aff_list *list2)
1896 {
1897         return pw_aff_list_set(list1, list2, &isl_pw_aff_ne_set);
1898 }
1899
1900 /* Return a set containing those elements in the shared domain
1901  * of the elements of list1 and list2 where each element in list1
1902  * is less than or equal to each element in list2.
1903  */
1904 __isl_give isl_set *isl_pw_aff_list_le_set(__isl_take isl_pw_aff_list *list1,
1905         __isl_take isl_pw_aff_list *list2)
1906 {
1907         return pw_aff_list_set(list1, list2, &isl_pw_aff_le_set);
1908 }
1909
1910 __isl_give isl_set *isl_pw_aff_list_lt_set(__isl_take isl_pw_aff_list *list1,
1911         __isl_take isl_pw_aff_list *list2)
1912 {
1913         return pw_aff_list_set(list1, list2, &isl_pw_aff_lt_set);
1914 }
1915
1916 __isl_give isl_set *isl_pw_aff_list_ge_set(__isl_take isl_pw_aff_list *list1,
1917         __isl_take isl_pw_aff_list *list2)
1918 {
1919         return pw_aff_list_set(list1, list2, &isl_pw_aff_ge_set);
1920 }
1921
1922 __isl_give isl_set *isl_pw_aff_list_gt_set(__isl_take isl_pw_aff_list *list1,
1923         __isl_take isl_pw_aff_list *list2)
1924 {
1925         return pw_aff_list_set(list1, list2, &isl_pw_aff_gt_set);
1926 }
1927
1928
1929 /* Return a set containing those elements in the shared domain
1930  * of pwaff1 and pwaff2 where pwaff1 is not equal to pwaff2.
1931  */
1932 static __isl_give isl_set *pw_aff_ne_set(__isl_take isl_pw_aff *pwaff1,
1933         __isl_take isl_pw_aff *pwaff2)
1934 {
1935         isl_set *set_lt, *set_gt;
1936
1937         set_lt = isl_pw_aff_lt_set(isl_pw_aff_copy(pwaff1),
1938                                    isl_pw_aff_copy(pwaff2));
1939         set_gt = isl_pw_aff_gt_set(pwaff1, pwaff2);
1940         return isl_set_union_disjoint(set_lt, set_gt);
1941 }
1942
1943 __isl_give isl_set *isl_pw_aff_ne_set(__isl_take isl_pw_aff *pwaff1,
1944         __isl_take isl_pw_aff *pwaff2)
1945 {
1946         return align_params_pw_pw_set_and(pwaff1, pwaff2, &pw_aff_ne_set);
1947 }
1948
1949 __isl_give isl_pw_aff *isl_pw_aff_scale_down(__isl_take isl_pw_aff *pwaff,
1950         isl_int v)
1951 {
1952         int i;
1953
1954         if (isl_int_is_one(v))
1955                 return pwaff;
1956         if (!isl_int_is_pos(v))
1957                 isl_die(isl_pw_aff_get_ctx(pwaff), isl_error_invalid,
1958                         "factor needs to be positive",
1959                         return isl_pw_aff_free(pwaff));
1960         pwaff = isl_pw_aff_cow(pwaff);
1961         if (!pwaff)
1962                 return NULL;
1963         if (pwaff->n == 0)
1964                 return pwaff;
1965
1966         for (i = 0; i < pwaff->n; ++i) {
1967                 pwaff->p[i].aff = isl_aff_scale_down(pwaff->p[i].aff, v);
1968                 if (!pwaff->p[i].aff)
1969                         return isl_pw_aff_free(pwaff);
1970         }
1971
1972         return pwaff;
1973 }
1974
1975 __isl_give isl_pw_aff *isl_pw_aff_floor(__isl_take isl_pw_aff *pwaff)
1976 {
1977         int i;
1978
1979         pwaff = isl_pw_aff_cow(pwaff);
1980         if (!pwaff)
1981                 return NULL;
1982         if (pwaff->n == 0)
1983                 return pwaff;
1984
1985         for (i = 0; i < pwaff->n; ++i) {
1986                 pwaff->p[i].aff = isl_aff_floor(pwaff->p[i].aff);
1987                 if (!pwaff->p[i].aff)
1988                         return isl_pw_aff_free(pwaff);
1989         }
1990
1991         return pwaff;
1992 }
1993
1994 __isl_give isl_pw_aff *isl_pw_aff_ceil(__isl_take isl_pw_aff *pwaff)
1995 {
1996         int i;
1997
1998         pwaff = isl_pw_aff_cow(pwaff);
1999         if (!pwaff)
2000                 return NULL;
2001         if (pwaff->n == 0)
2002                 return pwaff;
2003
2004         for (i = 0; i < pwaff->n; ++i) {
2005                 pwaff->p[i].aff = isl_aff_ceil(pwaff->p[i].aff);
2006                 if (!pwaff->p[i].aff)
2007                         return isl_pw_aff_free(pwaff);
2008         }
2009
2010         return pwaff;
2011 }
2012
2013 /* Assuming that "cond1" and "cond2" are disjoint,
2014  * return an affine expression that is equal to pwaff1 on cond1
2015  * and to pwaff2 on cond2.
2016  */
2017 static __isl_give isl_pw_aff *isl_pw_aff_select(
2018         __isl_take isl_set *cond1, __isl_take isl_pw_aff *pwaff1,
2019         __isl_take isl_set *cond2, __isl_take isl_pw_aff *pwaff2)
2020 {
2021         pwaff1 = isl_pw_aff_intersect_domain(pwaff1, cond1);
2022         pwaff2 = isl_pw_aff_intersect_domain(pwaff2, cond2);
2023
2024         return isl_pw_aff_add_disjoint(pwaff1, pwaff2);
2025 }
2026
2027 /* Return an affine expression that is equal to pwaff_true for elements
2028  * where "cond" is non-zero and to pwaff_false for elements where "cond"
2029  * is zero.
2030  * That is, return cond ? pwaff_true : pwaff_false;
2031  */
2032 __isl_give isl_pw_aff *isl_pw_aff_cond(__isl_take isl_pw_aff *cond,
2033         __isl_take isl_pw_aff *pwaff_true, __isl_take isl_pw_aff *pwaff_false)
2034 {
2035         isl_set *cond_true, *cond_false;
2036
2037         cond_true = isl_pw_aff_non_zero_set(isl_pw_aff_copy(cond));
2038         cond_false = isl_pw_aff_zero_set(cond);
2039         return isl_pw_aff_select(cond_true, pwaff_true,
2040                                  cond_false, pwaff_false);
2041 }
2042
2043 int isl_aff_is_cst(__isl_keep isl_aff *aff)
2044 {
2045         if (!aff)
2046                 return -1;
2047
2048         return isl_seq_first_non_zero(aff->v->el + 2, aff->v->size - 2) == -1;
2049 }
2050
2051 /* Check whether pwaff is a piecewise constant.
2052  */
2053 int isl_pw_aff_is_cst(__isl_keep isl_pw_aff *pwaff)
2054 {
2055         int i;
2056
2057         if (!pwaff)
2058                 return -1;
2059
2060         for (i = 0; i < pwaff->n; ++i) {
2061                 int is_cst = isl_aff_is_cst(pwaff->p[i].aff);
2062                 if (is_cst < 0 || !is_cst)
2063                         return is_cst;
2064         }
2065
2066         return 1;
2067 }
2068
2069 __isl_give isl_aff *isl_aff_mul(__isl_take isl_aff *aff1,
2070         __isl_take isl_aff *aff2)
2071 {
2072         if (!isl_aff_is_cst(aff2) && isl_aff_is_cst(aff1))
2073                 return isl_aff_mul(aff2, aff1);
2074
2075         if (!isl_aff_is_cst(aff2))
2076                 isl_die(isl_aff_get_ctx(aff1), isl_error_invalid,
2077                         "at least one affine expression should be constant",
2078                         goto error);
2079
2080         aff1 = isl_aff_cow(aff1);
2081         if (!aff1 || !aff2)
2082                 goto error;
2083
2084         aff1 = isl_aff_scale(aff1, aff2->v->el[1]);
2085         aff1 = isl_aff_scale_down(aff1, aff2->v->el[0]);
2086
2087         isl_aff_free(aff2);
2088         return aff1;
2089 error:
2090         isl_aff_free(aff1);
2091         isl_aff_free(aff2);
2092         return NULL;
2093 }
2094
2095 static __isl_give isl_pw_aff *pw_aff_add(__isl_take isl_pw_aff *pwaff1,
2096         __isl_take isl_pw_aff *pwaff2)
2097 {
2098         return isl_pw_aff_on_shared_domain(pwaff1, pwaff2, &isl_aff_add);
2099 }
2100
2101 __isl_give isl_pw_aff *isl_pw_aff_add(__isl_take isl_pw_aff *pwaff1,
2102         __isl_take isl_pw_aff *pwaff2)
2103 {
2104         return isl_pw_aff_align_params_pw_pw_and(pwaff1, pwaff2, &pw_aff_add);
2105 }
2106
2107 __isl_give isl_pw_aff *isl_pw_aff_union_add(__isl_take isl_pw_aff *pwaff1,
2108         __isl_take isl_pw_aff *pwaff2)
2109 {
2110         return isl_pw_aff_union_add_(pwaff1, pwaff2);
2111 }
2112
2113 static __isl_give isl_pw_aff *pw_aff_mul(__isl_take isl_pw_aff *pwaff1,
2114         __isl_take isl_pw_aff *pwaff2)
2115 {
2116         return isl_pw_aff_on_shared_domain(pwaff1, pwaff2, &isl_aff_mul);
2117 }
2118
2119 __isl_give isl_pw_aff *isl_pw_aff_mul(__isl_take isl_pw_aff *pwaff1,
2120         __isl_take isl_pw_aff *pwaff2)
2121 {
2122         return isl_pw_aff_align_params_pw_pw_and(pwaff1, pwaff2, &pw_aff_mul);
2123 }
2124
2125 static __isl_give isl_pw_aff *pw_aff_min(__isl_take isl_pw_aff *pwaff1,
2126         __isl_take isl_pw_aff *pwaff2)
2127 {
2128         isl_set *le;
2129         isl_set *dom;
2130
2131         dom = isl_set_intersect(isl_pw_aff_domain(isl_pw_aff_copy(pwaff1)),
2132                                 isl_pw_aff_domain(isl_pw_aff_copy(pwaff2)));
2133         le = isl_pw_aff_le_set(isl_pw_aff_copy(pwaff1),
2134                                 isl_pw_aff_copy(pwaff2));
2135         dom = isl_set_subtract(dom, isl_set_copy(le));
2136         return isl_pw_aff_select(le, pwaff1, dom, pwaff2);
2137 }
2138
2139 __isl_give isl_pw_aff *isl_pw_aff_min(__isl_take isl_pw_aff *pwaff1,
2140         __isl_take isl_pw_aff *pwaff2)
2141 {
2142         return isl_pw_aff_align_params_pw_pw_and(pwaff1, pwaff2, &pw_aff_min);
2143 }
2144
2145 static __isl_give isl_pw_aff *pw_aff_max(__isl_take isl_pw_aff *pwaff1,
2146         __isl_take isl_pw_aff *pwaff2)
2147 {
2148         isl_set *ge;
2149         isl_set *dom;
2150
2151         dom = isl_set_intersect(isl_pw_aff_domain(isl_pw_aff_copy(pwaff1)),
2152                                 isl_pw_aff_domain(isl_pw_aff_copy(pwaff2)));
2153         ge = isl_pw_aff_ge_set(isl_pw_aff_copy(pwaff1),
2154                                 isl_pw_aff_copy(pwaff2));
2155         dom = isl_set_subtract(dom, isl_set_copy(ge));
2156         return isl_pw_aff_select(ge, pwaff1, dom, pwaff2);
2157 }
2158
2159 __isl_give isl_pw_aff *isl_pw_aff_max(__isl_take isl_pw_aff *pwaff1,
2160         __isl_take isl_pw_aff *pwaff2)
2161 {
2162         return isl_pw_aff_align_params_pw_pw_and(pwaff1, pwaff2, &pw_aff_max);
2163 }
2164
2165 static __isl_give isl_pw_aff *pw_aff_list_reduce(
2166         __isl_take isl_pw_aff_list *list,
2167         __isl_give isl_pw_aff *(*fn)(__isl_take isl_pw_aff *pwaff1,
2168                                         __isl_take isl_pw_aff *pwaff2))
2169 {
2170         int i;
2171         isl_ctx *ctx;
2172         isl_pw_aff *res;
2173
2174         if (!list)
2175                 return NULL;
2176
2177         ctx = isl_pw_aff_list_get_ctx(list);
2178         if (list->n < 1)
2179                 isl_die(ctx, isl_error_invalid,
2180                         "list should contain at least one element",
2181                         return isl_pw_aff_list_free(list));
2182
2183         res = isl_pw_aff_copy(list->p[0]);
2184         for (i = 1; i < list->n; ++i)
2185                 res = fn(res, isl_pw_aff_copy(list->p[i]));
2186
2187         isl_pw_aff_list_free(list);
2188         return res;
2189 }
2190
2191 /* Return an isl_pw_aff that maps each element in the intersection of the
2192  * domains of the elements of list to the minimal corresponding affine
2193  * expression.
2194  */
2195 __isl_give isl_pw_aff *isl_pw_aff_list_min(__isl_take isl_pw_aff_list *list)
2196 {
2197         return pw_aff_list_reduce(list, &isl_pw_aff_min);
2198 }
2199
2200 /* Return an isl_pw_aff that maps each element in the intersection of the
2201  * domains of the elements of list to the maximal corresponding affine
2202  * expression.
2203  */
2204 __isl_give isl_pw_aff *isl_pw_aff_list_max(__isl_take isl_pw_aff_list *list)
2205 {
2206         return pw_aff_list_reduce(list, &isl_pw_aff_max);
2207 }
2208
2209 #undef BASE
2210 #define BASE aff
2211
2212 #include <isl_multi_templ.c>
2213
2214 /* Construct an isl_multi_aff in the given space with value zero in
2215  * each of the output dimensions.
2216  */
2217 __isl_give isl_multi_aff *isl_multi_aff_zero(__isl_take isl_space *space)
2218 {
2219         int n;
2220         isl_multi_aff *ma;
2221
2222         if (!space)
2223                 return NULL;
2224
2225         n = isl_space_dim(space , isl_dim_out);
2226         ma = isl_multi_aff_alloc(isl_space_copy(space));
2227
2228         if (!n)
2229                 isl_space_free(space);
2230         else {
2231                 int i;
2232                 isl_local_space *ls;
2233                 isl_aff *aff;
2234
2235                 space = isl_space_domain(space);
2236                 ls = isl_local_space_from_space(space);
2237                 aff = isl_aff_zero_on_domain(ls);
2238
2239                 for (i = 0; i < n; ++i)
2240                         ma = isl_multi_aff_set_aff(ma, i, isl_aff_copy(aff));
2241
2242                 isl_aff_free(aff);
2243         }
2244
2245         return ma;
2246 }
2247
2248 /* Create an isl_multi_aff in the given space that maps each
2249  * input dimension to the corresponding output dimension.
2250  */
2251 __isl_give isl_multi_aff *isl_multi_aff_identity(__isl_take isl_space *space)
2252 {
2253         int n;
2254         isl_multi_aff *ma;
2255
2256         if (!space)
2257                 return NULL;
2258
2259         if (isl_space_is_set(space))
2260                 isl_die(isl_space_get_ctx(space), isl_error_invalid,
2261                         "expecting map space", goto error);
2262
2263         n = isl_space_dim(space, isl_dim_out);
2264         if (n != isl_space_dim(space, isl_dim_in))
2265                 isl_die(isl_space_get_ctx(space), isl_error_invalid,
2266                         "number of input and output dimensions needs to be "
2267                         "the same", goto error);
2268
2269         ma = isl_multi_aff_alloc(isl_space_copy(space));
2270
2271         if (!n)
2272                 isl_space_free(space);
2273         else {
2274                 int i;
2275                 isl_local_space *ls;
2276                 isl_aff *aff;
2277
2278                 space = isl_space_domain(space);
2279                 ls = isl_local_space_from_space(space);
2280                 aff = isl_aff_zero_on_domain(ls);
2281
2282                 for (i = 0; i < n; ++i) {
2283                         isl_aff *aff_i;
2284                         aff_i = isl_aff_copy(aff);
2285                         aff_i = isl_aff_add_coefficient_si(aff_i,
2286                                                             isl_dim_in, i, 1);
2287                         ma = isl_multi_aff_set_aff(ma, i, aff_i);
2288                 }
2289
2290                 isl_aff_free(aff);
2291         }
2292
2293         return ma;
2294 error:
2295         isl_space_free(space);
2296         return NULL;
2297 }
2298
2299 /* Create an isl_pw_multi_aff with the given isl_multi_aff on a universe
2300  * domain.
2301  */
2302 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_from_multi_aff(
2303         __isl_take isl_multi_aff *ma)
2304 {
2305         isl_set *dom = isl_set_universe(isl_multi_aff_get_domain_space(ma));
2306         return isl_pw_multi_aff_alloc(dom, ma);
2307 }
2308
2309 /* Create a piecewise multi-affine expression in the given space that maps each
2310  * input dimension to the corresponding output dimension.
2311  */
2312 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_identity(
2313         __isl_take isl_space *space)
2314 {
2315         return isl_pw_multi_aff_from_multi_aff(isl_multi_aff_identity(space));
2316 }
2317
2318 __isl_give isl_multi_aff *isl_multi_aff_add(__isl_take isl_multi_aff *maff1,
2319         __isl_take isl_multi_aff *maff2)
2320 {
2321         int i;
2322         isl_ctx *ctx;
2323
2324         maff1 = isl_multi_aff_cow(maff1);
2325         if (!maff1 || !maff2)
2326                 goto error;
2327
2328         ctx = isl_multi_aff_get_ctx(maff1);
2329         if (!isl_space_is_equal(maff1->space, maff2->space))
2330                 isl_die(ctx, isl_error_invalid,
2331                         "spaces don't match", goto error);
2332
2333         for (i = 0; i < maff1->n; ++i) {
2334                 maff1->p[i] = isl_aff_add(maff1->p[i],
2335                                             isl_aff_copy(maff2->p[i]));
2336                 if (!maff1->p[i])
2337                         goto error;
2338         }
2339
2340         isl_multi_aff_free(maff2);
2341         return maff1;
2342 error:
2343         isl_multi_aff_free(maff1);
2344         isl_multi_aff_free(maff2);
2345         return NULL;
2346 }
2347
2348 /* Given two multi-affine expressions A -> B and C -> D,
2349  * construct a multi-affine expression [A -> C] -> [B -> D].
2350  */
2351 __isl_give isl_multi_aff *isl_multi_aff_product(
2352         __isl_take isl_multi_aff *ma1, __isl_take isl_multi_aff *ma2)
2353 {
2354         int i;
2355         isl_aff *aff;
2356         isl_space *space;
2357         isl_multi_aff *res;
2358         int in1, in2, out1, out2;
2359
2360         in1 = isl_multi_aff_dim(ma1, isl_dim_in);
2361         in2 = isl_multi_aff_dim(ma2, isl_dim_in);
2362         out1 = isl_multi_aff_dim(ma1, isl_dim_out);
2363         out2 = isl_multi_aff_dim(ma2, isl_dim_out);
2364         space = isl_space_product(isl_multi_aff_get_space(ma1),
2365                                   isl_multi_aff_get_space(ma2));
2366         res = isl_multi_aff_alloc(isl_space_copy(space));
2367         space = isl_space_domain(space);
2368
2369         for (i = 0; i < out1; ++i) {
2370                 aff = isl_multi_aff_get_aff(ma1, i);
2371                 aff = isl_aff_insert_dims(aff, isl_dim_in, in1, in2);
2372                 aff = isl_aff_reset_domain_space(aff, isl_space_copy(space));
2373                 res = isl_multi_aff_set_aff(res, i, aff);
2374         }
2375
2376         for (i = 0; i < out2; ++i) {
2377                 aff = isl_multi_aff_get_aff(ma2, i);
2378                 aff = isl_aff_insert_dims(aff, isl_dim_in, 0, in1);
2379                 aff = isl_aff_reset_domain_space(aff, isl_space_copy(space));
2380                 res = isl_multi_aff_set_aff(res, out1 + i, aff);
2381         }
2382
2383         isl_space_free(space);
2384         isl_multi_aff_free(ma1);
2385         isl_multi_aff_free(ma2);
2386         return res;
2387 }
2388
2389 /* Exploit the equalities in "eq" to simplify the affine expressions.
2390  */
2391 static __isl_give isl_multi_aff *isl_multi_aff_substitute_equalities(
2392         __isl_take isl_multi_aff *maff, __isl_take isl_basic_set *eq)
2393 {
2394         int i;
2395
2396         maff = isl_multi_aff_cow(maff);
2397         if (!maff || !eq)
2398                 goto error;
2399
2400         for (i = 0; i < maff->n; ++i) {
2401                 maff->p[i] = isl_aff_substitute_equalities(maff->p[i],
2402                                                     isl_basic_set_copy(eq));
2403                 if (!maff->p[i])
2404                         goto error;
2405         }
2406
2407         isl_basic_set_free(eq);
2408         return maff;
2409 error:
2410         isl_basic_set_free(eq);
2411         isl_multi_aff_free(maff);
2412         return NULL;
2413 }
2414
2415 __isl_give isl_multi_aff *isl_multi_aff_scale(__isl_take isl_multi_aff *maff,
2416         isl_int f)
2417 {
2418         int i;
2419
2420         maff = isl_multi_aff_cow(maff);
2421         if (!maff)
2422                 return NULL;
2423
2424         for (i = 0; i < maff->n; ++i) {
2425                 maff->p[i] = isl_aff_scale(maff->p[i], f);
2426                 if (!maff->p[i])
2427                         return isl_multi_aff_free(maff);
2428         }
2429
2430         return maff;
2431 }
2432
2433 __isl_give isl_multi_aff *isl_multi_aff_add_on_domain(__isl_keep isl_set *dom,
2434         __isl_take isl_multi_aff *maff1, __isl_take isl_multi_aff *maff2)
2435 {
2436         maff1 = isl_multi_aff_add(maff1, maff2);
2437         maff1 = isl_multi_aff_gist(maff1, isl_set_copy(dom));
2438         return maff1;
2439 }
2440
2441 int isl_multi_aff_is_empty(__isl_keep isl_multi_aff *maff)
2442 {
2443         if (!maff)
2444                 return -1;
2445
2446         return 0;
2447 }
2448
2449 int isl_multi_aff_plain_is_equal(__isl_keep isl_multi_aff *maff1,
2450         __isl_keep isl_multi_aff *maff2)
2451 {
2452         int i;
2453         int equal;
2454
2455         if (!maff1 || !maff2)
2456                 return -1;
2457         if (maff1->n != maff2->n)
2458                 return 0;
2459         equal = isl_space_is_equal(maff1->space, maff2->space);
2460         if (equal < 0 || !equal)
2461                 return equal;
2462
2463         for (i = 0; i < maff1->n; ++i) {
2464                 equal = isl_aff_plain_is_equal(maff1->p[i], maff2->p[i]);
2465                 if (equal < 0 || !equal)
2466                         return equal;
2467         }
2468
2469         return 1;
2470 }
2471
2472 __isl_give isl_multi_aff *isl_multi_aff_set_dim_name(
2473         __isl_take isl_multi_aff *maff,
2474         enum isl_dim_type type, unsigned pos, const char *s)
2475 {
2476         int i;
2477
2478         maff = isl_multi_aff_cow(maff);
2479         if (!maff)
2480                 return NULL;
2481
2482         maff->space = isl_space_set_dim_name(maff->space, type, pos, s);
2483         if (!maff->space)
2484                 return isl_multi_aff_free(maff);
2485
2486         if (type == isl_dim_out)
2487                 return maff;
2488         for (i = 0; i < maff->n; ++i) {
2489                 maff->p[i] = isl_aff_set_dim_name(maff->p[i], type, pos, s);
2490                 if (!maff->p[i])
2491                         return isl_multi_aff_free(maff);
2492         }
2493
2494         return maff;
2495 }
2496
2497 __isl_give isl_multi_aff *isl_multi_aff_drop_dims(__isl_take isl_multi_aff *maff,
2498         enum isl_dim_type type, unsigned first, unsigned n)
2499 {
2500         int i;
2501
2502         maff = isl_multi_aff_cow(maff);
2503         if (!maff)
2504                 return NULL;
2505
2506         maff->space = isl_space_drop_dims(maff->space, type, first, n);
2507         if (!maff->space)
2508                 return isl_multi_aff_free(maff);
2509
2510         if (type == isl_dim_out) {
2511                 for (i = 0; i < n; ++i)
2512                         isl_aff_free(maff->p[first + i]);
2513                 for (i = first; i + n < maff->n; ++i)
2514                         maff->p[i] = maff->p[i + n];
2515                 maff->n -= n;
2516                 return maff;
2517         }
2518
2519         for (i = 0; i < maff->n; ++i) {
2520                 maff->p[i] = isl_aff_drop_dims(maff->p[i], type, first, n);
2521                 if (!maff->p[i])
2522                         return isl_multi_aff_free(maff);
2523         }
2524
2525         return maff;
2526 }
2527
2528 /* Return the set of domain elements where "ma1" is lexicographically
2529  * smaller than or equal to "ma2".
2530  */
2531 __isl_give isl_set *isl_multi_aff_lex_le_set(__isl_take isl_multi_aff *ma1,
2532         __isl_take isl_multi_aff *ma2)
2533 {
2534         return isl_multi_aff_lex_ge_set(ma2, ma1);
2535 }
2536
2537 /* Return the set of domain elements where "ma1" is lexicographically
2538  * greater than or equal to "ma2".
2539  */
2540 __isl_give isl_set *isl_multi_aff_lex_ge_set(__isl_take isl_multi_aff *ma1,
2541         __isl_take isl_multi_aff *ma2)
2542 {
2543         isl_space *space;
2544         isl_map *map1, *map2;
2545         isl_map *map, *ge;
2546
2547         map1 = isl_map_from_multi_aff(ma1);
2548         map2 = isl_map_from_multi_aff(ma2);
2549         map = isl_map_range_product(map1, map2);
2550         space = isl_space_range(isl_map_get_space(map));
2551         space = isl_space_domain(isl_space_unwrap(space));
2552         ge = isl_map_lex_ge(space);
2553         map = isl_map_intersect_range(map, isl_map_wrap(ge));
2554
2555         return isl_map_domain(map);
2556 }
2557
2558 #undef PW
2559 #define PW isl_pw_multi_aff
2560 #undef EL
2561 #define EL isl_multi_aff
2562 #undef EL_IS_ZERO
2563 #define EL_IS_ZERO is_empty
2564 #undef ZERO
2565 #define ZERO empty
2566 #undef IS_ZERO
2567 #define IS_ZERO is_empty
2568 #undef FIELD
2569 #define FIELD maff
2570 #undef DEFAULT_IS_ZERO
2571 #define DEFAULT_IS_ZERO 0
2572
2573 #define NO_NEG
2574 #define NO_EVAL
2575 #define NO_OPT
2576 #define NO_INVOLVES_DIMS
2577 #define NO_MOVE_DIMS
2578 #define NO_INSERT_DIMS
2579 #define NO_LIFT
2580 #define NO_MORPH
2581
2582 #include <isl_pw_templ.c>
2583
2584 #undef UNION
2585 #define UNION isl_union_pw_multi_aff
2586 #undef PART
2587 #define PART isl_pw_multi_aff
2588 #undef PARTS
2589 #define PARTS pw_multi_aff
2590 #define ALIGN_DOMAIN
2591
2592 #define NO_EVAL
2593
2594 #include <isl_union_templ.c>
2595
2596 /* Given a function "cmp" that returns the set of elements where
2597  * "ma1" is "better" than "ma2", return the intersection of this
2598  * set with "dom1" and "dom2".
2599  */
2600 static __isl_give isl_set *shared_and_better(__isl_keep isl_set *dom1,
2601         __isl_keep isl_set *dom2, __isl_keep isl_multi_aff *ma1,
2602         __isl_keep isl_multi_aff *ma2,
2603         __isl_give isl_set *(*cmp)(__isl_take isl_multi_aff *ma1,
2604                                     __isl_take isl_multi_aff *ma2))
2605 {
2606         isl_set *common;
2607         isl_set *better;
2608         int is_empty;
2609
2610         common = isl_set_intersect(isl_set_copy(dom1), isl_set_copy(dom2));
2611         is_empty = isl_set_plain_is_empty(common);
2612         if (is_empty >= 0 && is_empty)
2613                 return common;
2614         if (is_empty < 0)
2615                 return isl_set_free(common);
2616         better = cmp(isl_multi_aff_copy(ma1), isl_multi_aff_copy(ma2));
2617         better = isl_set_intersect(common, better);
2618
2619         return better;
2620 }
2621
2622 /* Given a function "cmp" that returns the set of elements where
2623  * "ma1" is "better" than "ma2", return a piecewise multi affine
2624  * expression defined on the union of the definition domains
2625  * of "pma1" and "pma2" that maps to the "best" of "pma1" and
2626  * "pma2" on each cell.  If only one of the two input functions
2627  * is defined on a given cell, then it is considered the best.
2628  */
2629 static __isl_give isl_pw_multi_aff *pw_multi_aff_union_opt(
2630         __isl_take isl_pw_multi_aff *pma1,
2631         __isl_take isl_pw_multi_aff *pma2,
2632         __isl_give isl_set *(*cmp)(__isl_take isl_multi_aff *ma1,
2633                                     __isl_take isl_multi_aff *ma2))
2634 {
2635         int i, j, n;
2636         isl_pw_multi_aff *res = NULL;
2637         isl_ctx *ctx;
2638         isl_set *set = NULL;
2639
2640         if (!pma1 || !pma2)
2641                 goto error;
2642
2643         ctx = isl_space_get_ctx(pma1->dim);
2644         if (!isl_space_is_equal(pma1->dim, pma2->dim))
2645                 isl_die(ctx, isl_error_invalid,
2646                         "arguments should live in the same space", goto error);
2647
2648         if (isl_pw_multi_aff_is_empty(pma1)) {
2649                 isl_pw_multi_aff_free(pma1);
2650                 return pma2;
2651         }
2652
2653         if (isl_pw_multi_aff_is_empty(pma2)) {
2654                 isl_pw_multi_aff_free(pma2);
2655                 return pma1;
2656         }
2657
2658         n = 2 * (pma1->n + 1) * (pma2->n + 1);
2659         res = isl_pw_multi_aff_alloc_size(isl_space_copy(pma1->dim), n);
2660
2661         for (i = 0; i < pma1->n; ++i) {
2662                 set = isl_set_copy(pma1->p[i].set);
2663                 for (j = 0; j < pma2->n; ++j) {
2664                         isl_set *better;
2665                         int is_empty;
2666
2667                         better = shared_and_better(pma2->p[j].set,
2668                                         pma1->p[i].set, pma2->p[j].maff,
2669                                         pma1->p[i].maff, cmp);
2670                         is_empty = isl_set_plain_is_empty(better);
2671                         if (is_empty < 0 || is_empty) {
2672                                 isl_set_free(better);
2673                                 if (is_empty < 0)
2674                                         goto error;
2675                                 continue;
2676                         }
2677                         set = isl_set_subtract(set, isl_set_copy(better));
2678
2679                         res = isl_pw_multi_aff_add_piece(res, better,
2680                                         isl_multi_aff_copy(pma2->p[j].maff));
2681                 }
2682                 res = isl_pw_multi_aff_add_piece(res, set,
2683                                         isl_multi_aff_copy(pma1->p[i].maff));
2684         }
2685
2686         for (j = 0; j < pma2->n; ++j) {
2687                 set = isl_set_copy(pma2->p[j].set);
2688                 for (i = 0; i < pma1->n; ++i)
2689                         set = isl_set_subtract(set,
2690                                         isl_set_copy(pma1->p[i].set));
2691                 res = isl_pw_multi_aff_add_piece(res, set,
2692                                         isl_multi_aff_copy(pma2->p[j].maff));
2693         }
2694
2695         isl_pw_multi_aff_free(pma1);
2696         isl_pw_multi_aff_free(pma2);
2697
2698         return res;
2699 error:
2700         isl_pw_multi_aff_free(pma1);
2701         isl_pw_multi_aff_free(pma2);
2702         isl_set_free(set);
2703         return isl_pw_multi_aff_free(res);
2704 }
2705
2706 static __isl_give isl_pw_multi_aff *pw_multi_aff_union_lexmax(
2707         __isl_take isl_pw_multi_aff *pma1,
2708         __isl_take isl_pw_multi_aff *pma2)
2709 {
2710         return pw_multi_aff_union_opt(pma1, pma2, &isl_multi_aff_lex_ge_set);
2711 }
2712
2713 /* Given two piecewise multi affine expressions, return a piecewise
2714  * multi-affine expression defined on the union of the definition domains
2715  * of the inputs that is equal to the lexicographic maximum of the two
2716  * inputs on each cell.  If only one of the two inputs is defined on
2717  * a given cell, then it is considered to be the maximum.
2718  */
2719 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_union_lexmax(
2720         __isl_take isl_pw_multi_aff *pma1,
2721         __isl_take isl_pw_multi_aff *pma2)
2722 {
2723         return isl_pw_multi_aff_align_params_pw_pw_and(pma1, pma2,
2724                                                     &pw_multi_aff_union_lexmax);
2725 }
2726
2727 static __isl_give isl_pw_multi_aff *pw_multi_aff_union_lexmin(
2728         __isl_take isl_pw_multi_aff *pma1,
2729         __isl_take isl_pw_multi_aff *pma2)
2730 {
2731         return pw_multi_aff_union_opt(pma1, pma2, &isl_multi_aff_lex_le_set);
2732 }
2733
2734 /* Given two piecewise multi affine expressions, return a piecewise
2735  * multi-affine expression defined on the union of the definition domains
2736  * of the inputs that is equal to the lexicographic minimum of the two
2737  * inputs on each cell.  If only one of the two inputs is defined on
2738  * a given cell, then it is considered to be the minimum.
2739  */
2740 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_union_lexmin(
2741         __isl_take isl_pw_multi_aff *pma1,
2742         __isl_take isl_pw_multi_aff *pma2)
2743 {
2744         return isl_pw_multi_aff_align_params_pw_pw_and(pma1, pma2,
2745                                                     &pw_multi_aff_union_lexmin);
2746 }
2747
2748 static __isl_give isl_pw_multi_aff *pw_multi_aff_add(
2749         __isl_take isl_pw_multi_aff *pma1, __isl_take isl_pw_multi_aff *pma2)
2750 {
2751         return isl_pw_multi_aff_on_shared_domain(pma1, pma2,
2752                                                 &isl_multi_aff_add);
2753 }
2754
2755 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_add(
2756         __isl_take isl_pw_multi_aff *pma1, __isl_take isl_pw_multi_aff *pma2)
2757 {
2758         return isl_pw_multi_aff_align_params_pw_pw_and(pma1, pma2,
2759                                                 &pw_multi_aff_add);
2760 }
2761
2762 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_union_add(
2763         __isl_take isl_pw_multi_aff *pma1, __isl_take isl_pw_multi_aff *pma2)
2764 {
2765         return isl_pw_multi_aff_union_add_(pma1, pma2);
2766 }
2767
2768 /* Given two piecewise multi-affine expressions A -> B and C -> D,
2769  * construct a piecewise multi-affine expression [A -> C] -> [B -> D].
2770  */
2771 static __isl_give isl_pw_multi_aff *pw_multi_aff_product(
2772         __isl_take isl_pw_multi_aff *pma1, __isl_take isl_pw_multi_aff *pma2)
2773 {
2774         int i, j, n;
2775         isl_space *space;
2776         isl_pw_multi_aff *res;
2777
2778         if (!pma1 || !pma2)
2779                 goto error;
2780
2781         n = pma1->n * pma2->n;
2782         space = isl_space_product(isl_space_copy(pma1->dim),
2783                                   isl_space_copy(pma2->dim));
2784         res = isl_pw_multi_aff_alloc_size(space, n);
2785
2786         for (i = 0; i < pma1->n; ++i) {
2787                 for (j = 0; j < pma2->n; ++j) {
2788                         isl_set *domain;
2789                         isl_multi_aff *ma;
2790
2791                         domain = isl_set_product(isl_set_copy(pma1->p[i].set),
2792                                                  isl_set_copy(pma2->p[j].set));
2793                         ma = isl_multi_aff_product(
2794                                         isl_multi_aff_copy(pma1->p[i].maff),
2795                                         isl_multi_aff_copy(pma2->p[i].maff));
2796                         res = isl_pw_multi_aff_add_piece(res, domain, ma);
2797                 }
2798         }
2799
2800         isl_pw_multi_aff_free(pma1);
2801         isl_pw_multi_aff_free(pma2);
2802         return res;
2803 error:
2804         isl_pw_multi_aff_free(pma1);
2805         isl_pw_multi_aff_free(pma2);
2806         return NULL;
2807 }
2808
2809 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_product(
2810         __isl_take isl_pw_multi_aff *pma1, __isl_take isl_pw_multi_aff *pma2)
2811 {
2812         return isl_pw_multi_aff_align_params_pw_pw_and(pma1, pma2,
2813                                                 &pw_multi_aff_product);
2814 }
2815
2816 /* Construct a map mapping the domain of the piecewise multi-affine expression
2817  * to its range, with each dimension in the range equated to the
2818  * corresponding affine expression on its cell.
2819  */
2820 __isl_give isl_map *isl_map_from_pw_multi_aff(__isl_take isl_pw_multi_aff *pma)
2821 {
2822         int i;
2823         isl_map *map;
2824
2825         if (!pma)
2826                 return NULL;
2827
2828         map = isl_map_empty(isl_pw_multi_aff_get_space(pma));
2829
2830         for (i = 0; i < pma->n; ++i) {
2831                 isl_multi_aff *maff;
2832                 isl_basic_map *bmap;
2833                 isl_map *map_i;
2834
2835                 maff = isl_multi_aff_copy(pma->p[i].maff);
2836                 bmap = isl_basic_map_from_multi_aff(maff);
2837                 map_i = isl_map_from_basic_map(bmap);
2838                 map_i = isl_map_intersect_domain(map_i,
2839                                                 isl_set_copy(pma->p[i].set));
2840                 map = isl_map_union_disjoint(map, map_i);
2841         }
2842
2843         isl_pw_multi_aff_free(pma);
2844         return map;
2845 }
2846
2847 __isl_give isl_set *isl_set_from_pw_multi_aff(__isl_take isl_pw_multi_aff *pma)
2848 {
2849         if (!isl_space_is_set(pma->dim))
2850                 isl_die(isl_pw_multi_aff_get_ctx(pma), isl_error_invalid,
2851                         "isl_pw_multi_aff cannot be converted into an isl_set",
2852                         return isl_pw_multi_aff_free(pma));
2853
2854         return isl_map_from_pw_multi_aff(pma);
2855 }
2856
2857 /* Given a basic map with a single output dimension that is defined
2858  * in terms of the parameters and input dimensions using an equality,
2859  * extract an isl_aff that expresses the output dimension in terms
2860  * of the parameters and input dimensions.
2861  *
2862  * Since some applications expect the result of isl_pw_multi_aff_from_map
2863  * to only contain integer affine expressions, we compute the floor
2864  * of the expression before returning.
2865  *
2866  * This function shares some similarities with
2867  * isl_basic_map_has_defining_equality and isl_constraint_get_bound.
2868  */
2869 static __isl_give isl_aff *extract_isl_aff_from_basic_map(
2870         __isl_take isl_basic_map *bmap)
2871 {
2872         int i;
2873         unsigned offset;
2874         unsigned total;
2875         isl_local_space *ls;
2876         isl_aff *aff;
2877
2878         if (!bmap)
2879                 return NULL;
2880         if (isl_basic_map_dim(bmap, isl_dim_out) != 1)
2881                 isl_die(isl_basic_map_get_ctx(bmap), isl_error_invalid,
2882                         "basic map should have a single output dimension",
2883                         goto error);
2884         offset = isl_basic_map_offset(bmap, isl_dim_out);
2885         total = isl_basic_map_total_dim(bmap);
2886         for (i = 0; i < bmap->n_eq; ++i) {
2887                 if (isl_int_is_zero(bmap->eq[i][offset]))
2888                         continue;
2889                 if (isl_seq_first_non_zero(bmap->eq[i] + offset + 1,
2890                                            1 + total - (offset + 1)) != -1)
2891                         continue;
2892                 break;
2893         }
2894         if (i >= bmap->n_eq)
2895                 isl_die(isl_basic_map_get_ctx(bmap), isl_error_invalid,
2896                         "unable to find suitable equality", goto error);
2897         ls = isl_basic_map_get_local_space(bmap);
2898         aff = isl_aff_alloc(isl_local_space_domain(ls));
2899         if (!aff)
2900                 goto error;
2901         if (isl_int_is_neg(bmap->eq[i][offset]))
2902                 isl_seq_cpy(aff->v->el + 1, bmap->eq[i], offset);
2903         else
2904                 isl_seq_neg(aff->v->el + 1, bmap->eq[i], offset);
2905         isl_seq_clr(aff->v->el + 1 + offset, aff->v->size - (1 + offset));
2906         isl_int_abs(aff->v->el[0], bmap->eq[i][offset]);
2907         isl_basic_map_free(bmap);
2908
2909         aff = isl_aff_remove_unused_divs(aff);
2910         aff = isl_aff_floor(aff);
2911         return aff;
2912 error:
2913         isl_basic_map_free(bmap);
2914         return NULL;
2915 }
2916
2917 /* Given a basic map where each output dimension is defined
2918  * in terms of the parameters and input dimensions using an equality,
2919  * extract an isl_multi_aff that expresses the output dimensions in terms
2920  * of the parameters and input dimensions.
2921  */
2922 static __isl_give isl_multi_aff *extract_isl_multi_aff_from_basic_map(
2923         __isl_take isl_basic_map *bmap)
2924 {
2925         int i;
2926         unsigned n_out;
2927         isl_multi_aff *ma;
2928
2929         if (!bmap)
2930                 return NULL;
2931
2932         ma = isl_multi_aff_alloc(isl_basic_map_get_space(bmap));
2933         n_out = isl_basic_map_dim(bmap, isl_dim_out);
2934
2935         for (i = 0; i < n_out; ++i) {
2936                 isl_basic_map *bmap_i;
2937                 isl_aff *aff;
2938
2939                 bmap_i = isl_basic_map_copy(bmap);
2940                 bmap_i = isl_basic_map_project_out(bmap_i, isl_dim_out,
2941                                                         i + 1, n_out - (1 + i));
2942                 bmap_i = isl_basic_map_project_out(bmap_i, isl_dim_out, 0, i);
2943                 aff = extract_isl_aff_from_basic_map(bmap_i);
2944                 ma = isl_multi_aff_set_aff(ma, i, aff);
2945         }
2946
2947         isl_basic_map_free(bmap);
2948
2949         return ma;
2950 }
2951
2952 /* Create an isl_pw_multi_aff that is equivalent to
2953  * isl_map_intersect_domain(isl_map_from_basic_map(bmap), domain).
2954  * The given basic map is such that each output dimension is defined
2955  * in terms of the parameters and input dimensions using an equality.
2956  */
2957 static __isl_give isl_pw_multi_aff *plain_pw_multi_aff_from_map(
2958         __isl_take isl_set *domain, __isl_take isl_basic_map *bmap)
2959 {
2960         isl_multi_aff *ma;
2961
2962         ma = extract_isl_multi_aff_from_basic_map(bmap);
2963         return isl_pw_multi_aff_alloc(domain, ma);
2964 }
2965
2966 /* Try and create an isl_pw_multi_aff that is equivalent to the given isl_map.
2967  * This obivously only works if the input "map" is single-valued.
2968  * If so, we compute the lexicographic minimum of the image in the form
2969  * of an isl_pw_multi_aff.  Since the image is unique, it is equal
2970  * to its lexicographic minimum.
2971  * If the input is not single-valued, we produce an error.
2972  *
2973  * As a special case, we first check if all output dimensions are uniquely
2974  * defined in terms of the parameters and input dimensions over the entire
2975  * domain.  If so, we extract the desired isl_pw_multi_aff directly
2976  * from the affine hull of "map" and its domain.
2977  */
2978 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_from_map(__isl_take isl_map *map)
2979 {
2980         int i;
2981         int sv;
2982         isl_pw_multi_aff *pma;
2983         isl_basic_map *hull;
2984
2985         if (!map)
2986                 return NULL;
2987
2988         hull = isl_map_affine_hull(isl_map_copy(map));
2989         sv = isl_basic_map_plain_is_single_valued(hull);
2990         if (sv >= 0 && sv)
2991                 return plain_pw_multi_aff_from_map(isl_map_domain(map), hull);
2992         isl_basic_map_free(hull);
2993         if (sv < 0)
2994                 goto error;
2995
2996         sv = isl_map_is_single_valued(map);
2997         if (sv < 0)
2998                 goto error;
2999         if (!sv)
3000                 isl_die(isl_map_get_ctx(map), isl_error_invalid,
3001                         "map is not single-valued", goto error);
3002         map = isl_map_make_disjoint(map);
3003         if (!map)
3004                 return NULL;
3005
3006         pma = isl_pw_multi_aff_empty(isl_map_get_space(map));
3007
3008         for (i = 0; i < map->n; ++i) {
3009                 isl_pw_multi_aff *pma_i;
3010                 isl_basic_map *bmap;
3011                 bmap = isl_basic_map_copy(map->p[i]);
3012                 pma_i = isl_basic_map_lexmin_pw_multi_aff(bmap);
3013                 pma = isl_pw_multi_aff_add_disjoint(pma, pma_i);
3014         }
3015
3016         isl_map_free(map);
3017         return pma;
3018 error:
3019         isl_map_free(map);
3020         return NULL;
3021 }
3022
3023 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_from_set(__isl_take isl_set *set)
3024 {
3025         return isl_pw_multi_aff_from_map(set);
3026 }
3027
3028 /* Return the piecewise affine expression "set ? 1 : 0".
3029  */
3030 __isl_give isl_pw_aff *isl_set_indicator_function(__isl_take isl_set *set)
3031 {
3032         isl_pw_aff *pa;
3033         isl_space *space = isl_set_get_space(set);
3034         isl_local_space *ls = isl_local_space_from_space(space);
3035         isl_aff *zero = isl_aff_zero_on_domain(isl_local_space_copy(ls));
3036         isl_aff *one = isl_aff_zero_on_domain(ls);
3037
3038         one = isl_aff_add_constant_si(one, 1);
3039         pa = isl_pw_aff_alloc(isl_set_copy(set), one);
3040         set = isl_set_complement(set);
3041         pa = isl_pw_aff_add_disjoint(pa, isl_pw_aff_alloc(set, zero));
3042
3043         return pa;
3044 }
3045
3046 /* Plug in "subs" for dimension "type", "pos" of "aff".
3047  *
3048  * Let i be the dimension to replace and let "subs" be of the form
3049  *
3050  *      f/d
3051  *
3052  * and "aff" of the form
3053  *
3054  *      (a i + g)/m
3055  *
3056  * The result is
3057  *
3058  *      (a f + d g')/(m d)
3059  *
3060  * where g' is the result of plugging in "subs" in each of the integer
3061  * divisions in g.
3062  */
3063 __isl_give isl_aff *isl_aff_substitute(__isl_take isl_aff *aff,
3064         enum isl_dim_type type, unsigned pos, __isl_keep isl_aff *subs)
3065 {
3066         isl_ctx *ctx;
3067         isl_int v;
3068
3069         aff = isl_aff_cow(aff);
3070         if (!aff || !subs)
3071                 return isl_aff_free(aff);
3072
3073         ctx = isl_aff_get_ctx(aff);
3074         if (!isl_space_is_equal(aff->ls->dim, subs->ls->dim))
3075                 isl_die(ctx, isl_error_invalid,
3076                         "spaces don't match", return isl_aff_free(aff));
3077         if (isl_local_space_dim(subs->ls, isl_dim_div) != 0)
3078                 isl_die(ctx, isl_error_unsupported,
3079                         "cannot handle divs yet", return isl_aff_free(aff));
3080
3081         aff->ls = isl_local_space_substitute(aff->ls, type, pos, subs);
3082         if (!aff->ls)
3083                 return isl_aff_free(aff);
3084
3085         aff->v = isl_vec_cow(aff->v);
3086         if (!aff->v)
3087                 return isl_aff_free(aff);
3088
3089         pos += isl_local_space_offset(aff->ls, type);
3090
3091         isl_int_init(v);
3092         isl_int_set(v, aff->v->el[1 + pos]);
3093         isl_int_set_si(aff->v->el[1 + pos], 0);
3094         isl_seq_combine(aff->v->el + 1, subs->v->el[0], aff->v->el + 1,
3095                         v, subs->v->el + 1, subs->v->size - 1);
3096         isl_int_mul(aff->v->el[0], aff->v->el[0], subs->v->el[0]);
3097         isl_int_clear(v);
3098
3099         return aff;
3100 }
3101
3102 /* Plug in "subs" for dimension "type", "pos" in each of the affine
3103  * expressions in "maff".
3104  */
3105 __isl_give isl_multi_aff *isl_multi_aff_substitute(
3106         __isl_take isl_multi_aff *maff, enum isl_dim_type type, unsigned pos,
3107         __isl_keep isl_aff *subs)
3108 {
3109         int i;
3110
3111         maff = isl_multi_aff_cow(maff);
3112         if (!maff || !subs)
3113                 return isl_multi_aff_free(maff);
3114
3115         if (type == isl_dim_in)
3116                 type = isl_dim_set;
3117
3118         for (i = 0; i < maff->n; ++i) {
3119                 maff->p[i] = isl_aff_substitute(maff->p[i], type, pos, subs);
3120                 if (!maff->p[i])
3121                         return isl_multi_aff_free(maff);
3122         }
3123
3124         return maff;
3125 }
3126
3127 /* Plug in "subs" for dimension "type", "pos" of "pma".
3128  *
3129  * pma is of the form
3130  *
3131  *      A_i(v) -> M_i(v)
3132  *
3133  * while subs is of the form
3134  *
3135  *      v' = B_j(v) -> S_j
3136  *
3137  * Each pair i,j such that C_ij = A_i \cap B_i is non-empty
3138  * has a contribution in the result, in particular
3139  *
3140  *      C_ij(S_j) -> M_i(S_j)
3141  *
3142  * Note that plugging in S_j in C_ij may also result in an empty set
3143  * and this contribution should simply be discarded.
3144  */
3145 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_substitute(
3146         __isl_take isl_pw_multi_aff *pma, enum isl_dim_type type, unsigned pos,
3147         __isl_keep isl_pw_aff *subs)
3148 {
3149         int i, j, n;
3150         isl_pw_multi_aff *res;
3151
3152         if (!pma || !subs)
3153                 return isl_pw_multi_aff_free(pma);
3154
3155         n = pma->n * subs->n;
3156         res = isl_pw_multi_aff_alloc_size(isl_space_copy(pma->dim), n);
3157
3158         for (i = 0; i < pma->n; ++i) {
3159                 for (j = 0; j < subs->n; ++j) {
3160                         isl_set *common;
3161                         isl_multi_aff *res_ij;
3162                         common = isl_set_intersect(
3163                                         isl_set_copy(pma->p[i].set),
3164                                         isl_set_copy(subs->p[j].set));
3165                         common = isl_set_substitute(common,
3166                                         type, pos, subs->p[j].aff);
3167                         if (isl_set_plain_is_empty(common)) {
3168                                 isl_set_free(common);
3169                                 continue;
3170                         }
3171
3172                         res_ij = isl_multi_aff_substitute(
3173                                         isl_multi_aff_copy(pma->p[i].maff),
3174                                         type, pos, subs->p[j].aff);
3175
3176                         res = isl_pw_multi_aff_add_piece(res, common, res_ij);
3177                 }
3178         }
3179
3180         isl_pw_multi_aff_free(pma);
3181         return res;
3182 }
3183
3184 /* Extend the local space of "dst" to include the divs
3185  * in the local space of "src".
3186  */
3187 __isl_give isl_aff *isl_aff_align_divs(__isl_take isl_aff *dst,
3188         __isl_keep isl_aff *src)
3189 {
3190         isl_ctx *ctx;
3191         int *exp1 = NULL;
3192         int *exp2 = NULL;
3193         isl_mat *div;
3194
3195         if (!src || !dst)
3196                 return isl_aff_free(dst);
3197
3198         ctx = isl_aff_get_ctx(src);
3199         if (!isl_space_is_equal(src->ls->dim, dst->ls->dim))
3200                 isl_die(ctx, isl_error_invalid,
3201                         "spaces don't match", goto error);
3202
3203         if (src->ls->div->n_row == 0)
3204                 return dst;
3205
3206         exp1 = isl_alloc_array(ctx, int, src->ls->div->n_row);
3207         exp2 = isl_alloc_array(ctx, int, dst->ls->div->n_row);
3208         if (!exp1 || !exp2)
3209                 goto error;
3210
3211         div = isl_merge_divs(src->ls->div, dst->ls->div, exp1, exp2);
3212         dst = isl_aff_expand_divs(dst, div, exp2);
3213         free(exp1);
3214         free(exp2);
3215
3216         return dst;
3217 error:
3218         free(exp1);
3219         free(exp2);
3220         return isl_aff_free(dst);
3221 }
3222
3223 /* Adjust the local spaces of the affine expressions in "maff"
3224  * such that they all have the save divs.
3225  */
3226 __isl_give isl_multi_aff *isl_multi_aff_align_divs(
3227         __isl_take isl_multi_aff *maff)
3228 {
3229         int i;
3230
3231         if (!maff)
3232                 return NULL;
3233         if (maff->n == 0)
3234                 return maff;
3235         maff = isl_multi_aff_cow(maff);
3236         if (!maff)
3237                 return NULL;
3238
3239         for (i = 1; i < maff->n; ++i)
3240                 maff->p[0] = isl_aff_align_divs(maff->p[0], maff->p[i]);
3241         for (i = 1; i < maff->n; ++i) {
3242                 maff->p[i] = isl_aff_align_divs(maff->p[i], maff->p[0]);
3243                 if (!maff->p[i])
3244                         return isl_multi_aff_free(maff);
3245         }
3246
3247         return maff;
3248 }
3249
3250 __isl_give isl_aff *isl_aff_lift(__isl_take isl_aff *aff)
3251 {
3252         aff = isl_aff_cow(aff);
3253         if (!aff)
3254                 return NULL;
3255
3256         aff->ls = isl_local_space_lift(aff->ls);
3257         if (!aff->ls)
3258                 return isl_aff_free(aff);
3259
3260         return aff;
3261 }
3262
3263 /* Lift "maff" to a space with extra dimensions such that the result
3264  * has no more existentially quantified variables.
3265  * If "ls" is not NULL, then *ls is assigned the local space that lies
3266  * at the basis of the lifting applied to "maff".
3267  */
3268 __isl_give isl_multi_aff *isl_multi_aff_lift(__isl_take isl_multi_aff *maff,
3269         __isl_give isl_local_space **ls)
3270 {
3271         int i;
3272         isl_space *space;
3273         unsigned n_div;
3274
3275         if (ls)
3276                 *ls = NULL;
3277
3278         if (!maff)
3279                 return NULL;
3280
3281         if (maff->n == 0) {
3282                 if (ls) {
3283                         isl_space *space = isl_multi_aff_get_domain_space(maff);
3284                         *ls = isl_local_space_from_space(space);
3285                         if (!*ls)
3286                                 return isl_multi_aff_free(maff);
3287                 }
3288                 return maff;
3289         }
3290
3291         maff = isl_multi_aff_cow(maff);
3292         maff = isl_multi_aff_align_divs(maff);
3293         if (!maff)
3294                 return NULL;
3295
3296         n_div = isl_aff_dim(maff->p[0], isl_dim_div);
3297         space = isl_multi_aff_get_space(maff);
3298         space = isl_space_lift(isl_space_domain(space), n_div);
3299         space = isl_space_extend_domain_with_range(space,
3300                                                 isl_multi_aff_get_space(maff));
3301         if (!space)
3302                 return isl_multi_aff_free(maff);
3303         isl_space_free(maff->space);
3304         maff->space = space;
3305
3306         if (ls) {
3307                 *ls = isl_aff_get_domain_local_space(maff->p[0]);
3308                 if (!*ls)
3309                         return isl_multi_aff_free(maff);
3310         }
3311
3312         for (i = 0; i < maff->n; ++i) {
3313                 maff->p[i] = isl_aff_lift(maff->p[i]);
3314                 if (!maff->p[i])
3315                         goto error;
3316         }
3317
3318         return maff;
3319 error:
3320         if (ls)
3321                 isl_local_space_free(*ls);
3322         return isl_multi_aff_free(maff);
3323 }
3324
3325
3326 /* Extract an isl_pw_aff corresponding to output dimension "pos" of "pma".
3327  */
3328 __isl_give isl_pw_aff *isl_pw_multi_aff_get_pw_aff(
3329         __isl_keep isl_pw_multi_aff *pma, int pos)
3330 {
3331         int i;
3332         int n_out;
3333         isl_space *space;
3334         isl_pw_aff *pa;
3335
3336         if (!pma)
3337                 return NULL;
3338
3339         n_out = isl_pw_multi_aff_dim(pma, isl_dim_out);
3340         if (pos < 0 || pos >= n_out)
3341                 isl_die(isl_pw_multi_aff_get_ctx(pma), isl_error_invalid,
3342                         "index out of bounds", return NULL);
3343
3344         space = isl_pw_multi_aff_get_space(pma);
3345         space = isl_space_drop_dims(space, isl_dim_out,
3346                                     pos + 1, n_out - pos - 1);
3347         space = isl_space_drop_dims(space, isl_dim_out, 0, pos);
3348
3349         pa = isl_pw_aff_alloc_size(space, pma->n);
3350         for (i = 0; i < pma->n; ++i) {
3351                 isl_aff *aff;
3352                 aff = isl_multi_aff_get_aff(pma->p[i].maff, pos);
3353                 pa = isl_pw_aff_add_piece(pa, isl_set_copy(pma->p[i].set), aff);
3354         }
3355
3356         return pa;
3357 }
3358
3359 /* Return an isl_pw_multi_aff with the given "set" as domain and
3360  * an unnamed zero-dimensional range.
3361  */
3362 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_from_domain(
3363         __isl_take isl_set *set)
3364 {
3365         isl_multi_aff *ma;
3366         isl_space *space;
3367
3368         space = isl_set_get_space(set);
3369         space = isl_space_from_domain(space);
3370         ma = isl_multi_aff_zero(space);
3371         return isl_pw_multi_aff_alloc(set, ma);
3372 }
3373
3374 /* Add an isl_pw_multi_aff with the given "set" as domain and
3375  * an unnamed zero-dimensional range to *user.
3376  */
3377 static int add_pw_multi_aff_from_domain(__isl_take isl_set *set, void *user)
3378 {
3379         isl_union_pw_multi_aff **upma = user;
3380         isl_pw_multi_aff *pma;
3381
3382         pma = isl_pw_multi_aff_from_domain(set);
3383         *upma = isl_union_pw_multi_aff_add_pw_multi_aff(*upma, pma);
3384
3385         return 0;
3386 }
3387
3388 /* Return an isl_union_pw_multi_aff with the given "uset" as domain and
3389  * an unnamed zero-dimensional range.
3390  */
3391 __isl_give isl_union_pw_multi_aff *isl_union_pw_multi_aff_from_domain(
3392         __isl_take isl_union_set *uset)
3393 {
3394         isl_space *space;
3395         isl_union_pw_multi_aff *upma;
3396
3397         if (!uset)
3398                 return NULL;
3399
3400         space = isl_union_set_get_space(uset);
3401         upma = isl_union_pw_multi_aff_empty(space);
3402
3403         if (isl_union_set_foreach_set(uset,
3404                                     &add_pw_multi_aff_from_domain, &upma) < 0)
3405                 goto error;
3406
3407         isl_union_set_free(uset);
3408         return upma;
3409 error:
3410         isl_union_set_free(uset);
3411         isl_union_pw_multi_aff_free(upma);
3412         return NULL;
3413 }
3414
3415 /* Convert "pma" to an isl_map and add it to *umap.
3416  */
3417 static int map_from_pw_multi_aff(__isl_take isl_pw_multi_aff *pma, void *user)
3418 {
3419         isl_union_map **umap = user;
3420         isl_map *map;
3421
3422         map = isl_map_from_pw_multi_aff(pma);
3423         *umap = isl_union_map_add_map(*umap, map);
3424
3425         return 0;
3426 }
3427
3428 /* Construct a union map mapping the domain of the union
3429  * piecewise multi-affine expression to its range, with each dimension
3430  * in the range equated to the corresponding affine expression on its cell.
3431  */
3432 __isl_give isl_union_map *isl_union_map_from_union_pw_multi_aff(
3433         __isl_take isl_union_pw_multi_aff *upma)
3434 {
3435         isl_space *space;
3436         isl_union_map *umap;
3437
3438         if (!upma)
3439                 return NULL;
3440
3441         space = isl_union_pw_multi_aff_get_space(upma);
3442         umap = isl_union_map_empty(space);
3443
3444         if (isl_union_pw_multi_aff_foreach_pw_multi_aff(upma,
3445                                         &map_from_pw_multi_aff, &umap) < 0)
3446                 goto error;
3447
3448         isl_union_pw_multi_aff_free(upma);
3449         return umap;
3450 error:
3451         isl_union_pw_multi_aff_free(upma);
3452         isl_union_map_free(umap);
3453         return NULL;
3454 }
3455
3456 /* Local data for bin_entry and the callback "fn".
3457  */
3458 struct isl_union_pw_multi_aff_bin_data {
3459         isl_union_pw_multi_aff *upma2;
3460         isl_union_pw_multi_aff *res;
3461         isl_pw_multi_aff *pma;
3462         int (*fn)(void **entry, void *user);
3463 };
3464
3465 /* Given an isl_pw_multi_aff from upma1, store it in data->pma
3466  * and call data->fn for each isl_pw_multi_aff in data->upma2.
3467  */
3468 static int bin_entry(void **entry, void *user)
3469 {
3470         struct isl_union_pw_multi_aff_bin_data *data = user;
3471         isl_pw_multi_aff *pma = *entry;
3472
3473         data->pma = pma;
3474         if (isl_hash_table_foreach(data->upma2->dim->ctx, &data->upma2->table,
3475                                    data->fn, data) < 0)
3476                 return -1;
3477
3478         return 0;
3479 }
3480
3481 /* Call "fn" on each pair of isl_pw_multi_affs in "upma1" and "upma2".
3482  * The isl_pw_multi_aff from upma1 is stored in data->pma (where data is
3483  * passed as user field) and the isl_pw_multi_aff from upma2 is available
3484  * as *entry.  The callback should adjust data->res if desired.
3485  */
3486 static __isl_give isl_union_pw_multi_aff *bin_op(
3487         __isl_take isl_union_pw_multi_aff *upma1,
3488         __isl_take isl_union_pw_multi_aff *upma2,
3489         int (*fn)(void **entry, void *user))
3490 {
3491         isl_space *space;
3492         struct isl_union_pw_multi_aff_bin_data data = { NULL, NULL, NULL, fn };
3493
3494         space = isl_union_pw_multi_aff_get_space(upma2);
3495         upma1 = isl_union_pw_multi_aff_align_params(upma1, space);
3496         space = isl_union_pw_multi_aff_get_space(upma1);
3497         upma2 = isl_union_pw_multi_aff_align_params(upma2, space);
3498
3499         if (!upma1 || !upma2)
3500                 goto error;
3501
3502         data.upma2 = upma2;
3503         data.res = isl_union_pw_multi_aff_alloc(isl_space_copy(upma1->dim),
3504                                        upma1->table.n);
3505         if (isl_hash_table_foreach(upma1->dim->ctx, &upma1->table,
3506                                    &bin_entry, &data) < 0)
3507                 goto error;
3508
3509         isl_union_pw_multi_aff_free(upma1);
3510         isl_union_pw_multi_aff_free(upma2);
3511         return data.res;
3512 error:
3513         isl_union_pw_multi_aff_free(upma1);
3514         isl_union_pw_multi_aff_free(upma2);
3515         isl_union_pw_multi_aff_free(data.res);
3516         return NULL;
3517 }
3518
3519 /* Given two isl_multi_affs A -> B and C -> D,
3520  * construct an isl_multi_aff (A * C) -> (B, D).
3521  */
3522 __isl_give isl_multi_aff *isl_multi_aff_flat_range_product(
3523         __isl_take isl_multi_aff *ma1, __isl_take isl_multi_aff *ma2)
3524 {
3525         int i, n1, n2;
3526         isl_aff *aff;
3527         isl_space *space;
3528         isl_multi_aff *res;
3529
3530         if (!ma1 || !ma2)
3531                 goto error;
3532
3533         space = isl_space_range_product(isl_multi_aff_get_space(ma1),
3534                                         isl_multi_aff_get_space(ma2));
3535         space = isl_space_flatten_range(space);
3536         res = isl_multi_aff_alloc(space);
3537
3538         n1 = isl_multi_aff_dim(ma1, isl_dim_out);
3539         n2 = isl_multi_aff_dim(ma2, isl_dim_out);
3540
3541         for (i = 0; i < n1; ++i) {
3542                 aff = isl_multi_aff_get_aff(ma1, i);
3543                 res = isl_multi_aff_set_aff(res, i, aff);
3544         }
3545
3546         for (i = 0; i < n2; ++i) {
3547                 aff = isl_multi_aff_get_aff(ma2, i);
3548                 res = isl_multi_aff_set_aff(res, n1 + i, aff);
3549         }
3550
3551         isl_multi_aff_free(ma1);
3552         isl_multi_aff_free(ma2);
3553         return res;
3554 error:
3555         isl_multi_aff_free(ma1);
3556         isl_multi_aff_free(ma2);
3557         return NULL;
3558 }
3559
3560 /* Given two aligned isl_pw_multi_affs A -> B and C -> D,
3561  * construct an isl_pw_multi_aff (A * C) -> (B, D).
3562  */
3563 static __isl_give isl_pw_multi_aff *pw_multi_aff_flat_range_product(
3564         __isl_take isl_pw_multi_aff *pma1, __isl_take isl_pw_multi_aff *pma2)
3565 {
3566         isl_space *space;
3567
3568         space = isl_space_range_product(isl_pw_multi_aff_get_space(pma1),
3569                                         isl_pw_multi_aff_get_space(pma2));
3570         space = isl_space_flatten_range(space);
3571         return isl_pw_multi_aff_on_shared_domain_in(pma1, pma2, space,
3572                                             &isl_multi_aff_flat_range_product);
3573 }
3574
3575 /* Given two isl_pw_multi_affs A -> B and C -> D,
3576  * construct an isl_pw_multi_aff (A * C) -> (B, D).
3577  */
3578 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_flat_range_product(
3579         __isl_take isl_pw_multi_aff *pma1, __isl_take isl_pw_multi_aff *pma2)
3580 {
3581         return isl_pw_multi_aff_align_params_pw_pw_and(pma1, pma2,
3582                                             &pw_multi_aff_flat_range_product);
3583 }
3584
3585 /* If data->pma and *entry have the same domain space, then compute
3586  * their flat range product and the result to data->res.
3587  */
3588 static int flat_range_product_entry(void **entry, void *user)
3589 {
3590         struct isl_union_pw_multi_aff_bin_data *data = user;
3591         isl_pw_multi_aff *pma2 = *entry;
3592
3593         if (!isl_space_tuple_match(data->pma->dim, isl_dim_in,
3594                                  pma2->dim, isl_dim_in))
3595                 return 0;
3596
3597         pma2 = isl_pw_multi_aff_flat_range_product(
3598                                         isl_pw_multi_aff_copy(data->pma),
3599                                         isl_pw_multi_aff_copy(pma2));
3600
3601         data->res = isl_union_pw_multi_aff_add_pw_multi_aff(data->res, pma2);
3602
3603         return 0;
3604 }
3605
3606 /* Given two isl_union_pw_multi_affs A -> B and C -> D,
3607  * construct an isl_union_pw_multi_aff (A * C) -> (B, D).
3608  */
3609 __isl_give isl_union_pw_multi_aff *isl_union_pw_multi_aff_flat_range_product(
3610         __isl_take isl_union_pw_multi_aff *upma1,
3611         __isl_take isl_union_pw_multi_aff *upma2)
3612 {
3613         return bin_op(upma1, upma2, &flat_range_product_entry);
3614 }
3615
3616 /* Replace the affine expressions at position "pos" in "pma" by "pa".
3617  * The parameters are assumed to have been aligned.
3618  *
3619  * The implementation essentially performs an isl_pw_*_on_shared_domain,
3620  * except that it works on two different isl_pw_* types.
3621  */
3622 static __isl_give isl_pw_multi_aff *pw_multi_aff_set_pw_aff(
3623         __isl_take isl_pw_multi_aff *pma, unsigned pos,
3624         __isl_take isl_pw_aff *pa)
3625 {
3626         int i, j, n;
3627         isl_pw_multi_aff *res = NULL;
3628
3629         if (!pma || !pa)
3630                 goto error;
3631
3632         if (!isl_space_tuple_match(pma->dim, isl_dim_in, pa->dim, isl_dim_in))
3633                 isl_die(isl_pw_multi_aff_get_ctx(pma), isl_error_invalid,
3634                         "domains don't match", goto error);
3635         if (pos >= isl_pw_multi_aff_dim(pma, isl_dim_out))
3636                 isl_die(isl_pw_multi_aff_get_ctx(pma), isl_error_invalid,
3637                         "index out of bounds", goto error);
3638
3639         n = pma->n * pa->n;
3640         res = isl_pw_multi_aff_alloc_size(isl_pw_multi_aff_get_space(pma), n);
3641
3642         for (i = 0; i < pma->n; ++i) {
3643                 for (j = 0; j < pa->n; ++j) {
3644                         isl_set *common;
3645                         isl_multi_aff *res_ij;
3646                         int empty;
3647
3648                         common = isl_set_intersect(isl_set_copy(pma->p[i].set),
3649                                                    isl_set_copy(pa->p[j].set));
3650                         empty = isl_set_plain_is_empty(common);
3651                         if (empty < 0 || empty) {
3652                                 isl_set_free(common);
3653                                 if (empty < 0)
3654                                         goto error;
3655                                 continue;
3656                         }
3657
3658                         res_ij = isl_multi_aff_set_aff(
3659                                         isl_multi_aff_copy(pma->p[i].maff), pos,
3660                                         isl_aff_copy(pa->p[j].aff));
3661                         res_ij = isl_multi_aff_gist(res_ij,
3662                                         isl_set_copy(common));
3663
3664                         res = isl_pw_multi_aff_add_piece(res, common, res_ij);
3665                 }
3666         }
3667
3668         isl_pw_multi_aff_free(pma);
3669         isl_pw_aff_free(pa);
3670         return res;
3671 error:
3672         isl_pw_multi_aff_free(pma);
3673         isl_pw_aff_free(pa);
3674         return isl_pw_multi_aff_free(res);
3675 }
3676
3677 /* Replace the affine expressions at position "pos" in "pma" by "pa".
3678  */
3679 __isl_give isl_pw_multi_aff *isl_pw_multi_aff_set_pw_aff(
3680         __isl_take isl_pw_multi_aff *pma, unsigned pos,
3681         __isl_take isl_pw_aff *pa)
3682 {
3683         if (!pma || !pa)
3684                 goto error;
3685         if (isl_space_match(pma->dim, isl_dim_param, pa->dim, isl_dim_param))
3686                 return pw_multi_aff_set_pw_aff(pma, pos, pa);
3687         if (!isl_space_has_named_params(pma->dim) ||
3688             !isl_space_has_named_params(pa->dim))
3689                 isl_die(isl_pw_multi_aff_get_ctx(pma), isl_error_invalid,
3690                         "unaligned unnamed parameters", goto error);
3691         pma = isl_pw_multi_aff_align_params(pma, isl_pw_aff_get_space(pa));
3692         pa = isl_pw_aff_align_params(pa, isl_pw_multi_aff_get_space(pma));
3693         return pw_multi_aff_set_pw_aff(pma, pos, pa);
3694 error:
3695         isl_pw_multi_aff_free(pma);
3696         isl_pw_aff_free(pa);
3697         return NULL;
3698 }