isl_aff.c

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