isl_aff.c

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