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