isl_fold.c

   1 /*
   2  * Copyright 2010      INRIA Saclay
   3  *
   4  * Use of this software is governed by the MIT license
   5  *
   6  * Written by Sven Verdoolaege, INRIA Saclay - Ile-de-France,
   7  * Parc Club Orsay Universite, ZAC des vignes, 4 rue Jacques Monod,
   8  * 91893 Orsay, France
   9  */
  10
  11 #define ISL_DIM_H
  12 #include <isl_map_private.h>
  13 #include <isl_union_map_private.h>
  14 #include <isl_polynomial_private.h>
  15 #include <isl_point_private.h>
  16 #include <isl_space_private.h>
  17 #include <isl_lp_private.h>
  18 #include <isl_seq.h>
  19 #include <isl_mat_private.h>
  20 #include <isl_val_private.h>
  21 #include <isl_vec_private.h>
  22 #include <isl_config.h>
  23 #include <isl/deprecated/polynomial_int.h>
  24
  25 enum isl_fold isl_fold_type_negate(enum isl_fold type)
  26 {
  27         switch (type) {
  28         case isl_fold_min:
  29                 return isl_fold_max;
  30         case isl_fold_max:
  31                 return isl_fold_min;
  32         case isl_fold_list:
  33                 return isl_fold_list;
  34         }
  35
  36         isl_die(NULL, isl_error_internal, "unhandled isl_fold type", abort());
  37 }
  38
  39 static __isl_give isl_qpolynomial_fold *qpolynomial_fold_alloc(
  40         enum isl_fold type, __isl_take isl_space *dim, int n)
  41 {
  42         isl_qpolynomial_fold *fold;
  43
  44         if (!dim)
  45                 goto error;
  46
  47         isl_assert(dim->ctx, n >= 0, goto error);
  48         fold = isl_calloc(dim->ctx, struct isl_qpolynomial_fold,
  49                         sizeof(struct isl_qpolynomial_fold) +
  50                         (n - 1) * sizeof(struct isl_qpolynomial *));
  51         if (!fold)
  52                 goto error;
  53
  54         fold->ref = 1;
  55         fold->size = n;
  56         fold->n = 0;
  57         fold->type = type;
  58         fold->dim = dim;
  59
  60         return fold;
  61 error:
  62         isl_space_free(dim);
  63         return NULL;
  64 }
  65
  66 isl_ctx *isl_qpolynomial_fold_get_ctx(__isl_keep isl_qpolynomial_fold *fold)
  67 {
  68         return fold ? fold->dim->ctx : NULL;
  69 }
  70
  71 __isl_give isl_space *isl_qpolynomial_fold_get_domain_space(
  72         __isl_keep isl_qpolynomial_fold *fold)
  73 {
  74         return fold ? isl_space_copy(fold->dim) : NULL;
  75 }
  76
  77 __isl_give isl_space *isl_qpolynomial_fold_get_space(
  78         __isl_keep isl_qpolynomial_fold *fold)
  79 {
  80         isl_space *space;
  81         if (!fold)
  82                 return NULL;
  83         space = isl_space_copy(fold->dim);
  84         space = isl_space_from_domain(space);
  85         space = isl_space_add_dims(space, isl_dim_out, 1);
  86         return space;
  87 }
  88
  89 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_reset_domain_space(
  90         __isl_take isl_qpolynomial_fold *fold, __isl_take isl_space *dim)
  91 {
  92         int i;
  93
  94         fold = isl_qpolynomial_fold_cow(fold);
  95         if (!fold || !dim)
  96                 goto error;
  97
  98         for (i = 0; i < fold->n; ++i) {
  99                 fold->qp[i] = isl_qpolynomial_reset_domain_space(fold->qp[i],
 100                                                         isl_space_copy(dim));
 101                 if (!fold->qp[i])
 102                         goto error;
 103         }
 104
 105         isl_space_free(fold->dim);
 106         fold->dim = dim;
 107
 108         return fold;
 109 error:
 110         isl_qpolynomial_fold_free(fold);
 111         isl_space_free(dim);
 112         return NULL;
 113 }
 114
 115 /* Reset the space of "fold".  This function is called from isl_pw_templ.c
 116  * and doesn't know if the space of an element object is represented
 117  * directly or through its domain.  It therefore passes along both.
 118  */
 119 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_reset_space_and_domain(
 120         __isl_take isl_qpolynomial_fold *fold, __isl_take isl_space *space,
 121         __isl_take isl_space *domain)
 122 {
 123         isl_space_free(space);
 124         return isl_qpolynomial_fold_reset_domain_space(fold, domain);
 125 }
 126
 127 int isl_qpolynomial_fold_involves_dims(__isl_keep isl_qpolynomial_fold *fold,
 128         enum isl_dim_type type, unsigned first, unsigned n)
 129 {
 130         int i;
 131
 132         if (!fold)
 133                 return -1;
 134         if (fold->n == 0 || n == 0)
 135                 return 0;
 136
 137         for (i = 0; i < fold->n; ++i) {
 138                 int involves = isl_qpolynomial_involves_dims(fold->qp[i],
 139                                                             type, first, n);
 140                 if (involves < 0 || involves)
 141                         return involves;
 142         }
 143         return 0;
 144 }
 145
 146 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_set_dim_name(
 147         __isl_take isl_qpolynomial_fold *fold,
 148         enum isl_dim_type type, unsigned pos, const char *s)
 149 {
 150         int i;
 151
 152         fold = isl_qpolynomial_fold_cow(fold);
 153         if (!fold)
 154                 return NULL;
 155         fold->dim = isl_space_set_dim_name(fold->dim, type, pos, s);
 156         if (!fold->dim)
 157                 goto error;
 158
 159         for (i = 0; i < fold->n; ++i) {
 160                 fold->qp[i] = isl_qpolynomial_set_dim_name(fold->qp[i],
 161                                                             type, pos, s);
 162                 if (!fold->qp[i])
 163                         goto error;
 164         }
 165
 166         return fold;
 167 error:
 168         isl_qpolynomial_fold_free(fold);
 169         return NULL;
 170 }
 171
 172 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_drop_dims(
 173         __isl_take isl_qpolynomial_fold *fold,
 174         enum isl_dim_type type, unsigned first, unsigned n)
 175 {
 176         int i;
 177         enum isl_dim_type set_type;
 178
 179         if (!fold)
 180                 return NULL;
 181         if (n == 0)
 182                 return fold;
 183
 184         set_type = type == isl_dim_in ? isl_dim_set : type;
 185
 186         fold = isl_qpolynomial_fold_cow(fold);
 187         if (!fold)
 188                 return NULL;
 189         fold->dim = isl_space_drop_dims(fold->dim, set_type, first, n);
 190         if (!fold->dim)
 191                 goto error;
 192
 193         for (i = 0; i < fold->n; ++i) {
 194                 fold->qp[i] = isl_qpolynomial_drop_dims(fold->qp[i],
 195                                                             type, first, n);
 196                 if (!fold->qp[i])
 197                         goto error;
 198         }
 199
 200         return fold;
 201 error:
 202         isl_qpolynomial_fold_free(fold);
 203         return NULL;
 204 }
 205
 206 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_insert_dims(
 207         __isl_take isl_qpolynomial_fold *fold,
 208         enum isl_dim_type type, unsigned first, unsigned n)
 209 {
 210         int i;
 211
 212         if (!fold)
 213                 return NULL;
 214         if (n == 0 && !isl_space_is_named_or_nested(fold->dim, type))
 215                 return fold;
 216
 217         fold = isl_qpolynomial_fold_cow(fold);
 218         if (!fold)
 219                 return NULL;
 220         fold->dim = isl_space_insert_dims(fold->dim, type, first, n);
 221         if (!fold->dim)
 222                 goto error;
 223
 224         for (i = 0; i < fold->n; ++i) {
 225                 fold->qp[i] = isl_qpolynomial_insert_dims(fold->qp[i],
 226                                                             type, first, n);
 227                 if (!fold->qp[i])
 228                         goto error;
 229         }
 230
 231         return fold;
 232 error:
 233         isl_qpolynomial_fold_free(fold);
 234         return NULL;
 235 }
 236
 237 /* Determine the sign of the constant quasipolynomial "qp".
 238  *
 239  * Return
 240  *      -1 if qp <= 0
 241  *       1 if qp >= 0
 242  *       0 if unknown
 243  *
 244  * For qp == 0, we can return either -1 or 1.  In practice, we return 1.
 245  * For qp == NaN, the sign is undefined, so we return 0.
 246  */
 247 static int isl_qpolynomial_cst_sign(__isl_keep isl_qpolynomial *qp)
 248 {
 249         struct isl_upoly_cst *cst;
 250
 251         if (isl_qpolynomial_is_nan(qp))
 252                 return 0;
 253
 254         cst = isl_upoly_as_cst(qp->upoly);
 255         if (!cst)
 256                 return 0;
 257
 258         return isl_int_sgn(cst->n) < 0 ? -1 : 1;
 259 }
 260
 261 static int isl_qpolynomial_aff_sign(__isl_keep isl_set *set,
 262         __isl_keep isl_qpolynomial *qp)
 263 {
 264         enum isl_lp_result res;
 265         isl_vec *aff;
 266         isl_int opt;
 267         int sgn = 0;
 268
 269         aff = isl_qpolynomial_extract_affine(qp);
 270         if (!aff)
 271                 return 0;
 272
 273         isl_int_init(opt);
 274
 275         res = isl_set_solve_lp(set, 0, aff->el + 1, aff->el[0],
 276                                 &opt, NULL, NULL);
 277         if (res == isl_lp_error)
 278                 goto done;
 279         if (res == isl_lp_empty ||
 280             (res == isl_lp_ok && !isl_int_is_neg(opt))) {
 281                 sgn = 1;
 282                 goto done;
 283         }
 284
 285         res = isl_set_solve_lp(set, 1, aff->el + 1, aff->el[0],
 286                                 &opt, NULL, NULL);
 287         if (res == isl_lp_ok && !isl_int_is_pos(opt))
 288                 sgn = -1;
 289
 290 done:
 291         isl_int_clear(opt);
 292         isl_vec_free(aff);
 293         return sgn;
 294 }
 295
 296 /* Determine, if possible, the sign of the quasipolynomial "qp" on
 297  * the domain "set".
 298  *
 299  * If qp is a constant, then the problem is trivial.
 300  * If qp is linear, then we check if the minimum of the corresponding
 301  * affine constraint is non-negative or if the maximum is non-positive.
 302  *
 303  * Otherwise, we check if the outermost variable "v" has a lower bound "l"
 304  * in "set".  If so, we write qp(v,v') as
 305  *
 306  *      q(v,v') * (v - l) + r(v')
 307  *
 308  * if q(v,v') and r(v') have the same known sign, then the original
 309  * quasipolynomial has the same sign as well.
 310  *
 311  * Return
 312  *      -1 if qp <= 0
 313  *       1 if qp >= 0
 314  *       0 if unknown
 315  */
 316 static int isl_qpolynomial_sign(__isl_keep isl_set *set,
 317         __isl_keep isl_qpolynomial *qp)
 318 {
 319         int d;
 320         int i;
 321         int is;
 322         struct isl_upoly_rec *rec;
 323         isl_vec *v;
 324         isl_int l;
 325         enum isl_lp_result res;
 326         int sgn = 0;
 327
 328         is = isl_qpolynomial_is_cst(qp, NULL, NULL);
 329         if (is < 0)
 330                 return 0;
 331         if (is)
 332                 return isl_qpolynomial_cst_sign(qp);
 333
 334         is = isl_qpolynomial_is_affine(qp);
 335         if (is < 0)
 336                 return 0;
 337         if (is)
 338                 return isl_qpolynomial_aff_sign(set, qp);
 339
 340         if (qp->div->n_row > 0)
 341                 return 0;
 342
 343         rec = isl_upoly_as_rec(qp->upoly);
 344         if (!rec)
 345                 return 0;
 346
 347         d = isl_space_dim(qp->dim, isl_dim_all);
 348         v = isl_vec_alloc(set->ctx, 2 + d);
 349         if (!v)
 350                 return 0;
 351
 352         isl_seq_clr(v->el + 1, 1 + d);
 353         isl_int_set_si(v->el[0], 1);
 354         isl_int_set_si(v->el[2 + qp->upoly->var], 1);
 355
 356         isl_int_init(l);
 357
 358         res = isl_set_solve_lp(set, 0, v->el + 1, v->el[0], &l, NULL, NULL);
 359         if (res == isl_lp_ok) {
 360                 isl_qpolynomial *min;
 361                 isl_qpolynomial *base;
 362                 isl_qpolynomial *r, *q;
 363                 isl_qpolynomial *t;
 364
 365                 min = isl_qpolynomial_cst_on_domain(isl_space_copy(qp->dim), l);
 366                 base = isl_qpolynomial_var_pow_on_domain(isl_space_copy(qp->dim),
 367                                                 qp->upoly->var, 1);
 368
 369                 r = isl_qpolynomial_alloc(isl_space_copy(qp->dim), 0,
 370                                           isl_upoly_copy(rec->p[rec->n - 1]));
 371                 q = isl_qpolynomial_copy(r);
 372
 373                 for (i = rec->n - 2; i >= 0; --i) {
 374                         r = isl_qpolynomial_mul(r, isl_qpolynomial_copy(min));
 375                         t = isl_qpolynomial_alloc(isl_space_copy(qp->dim), 0,
 376                                                   isl_upoly_copy(rec->p[i]));
 377                         r = isl_qpolynomial_add(r, t);
 378                         if (i == 0)
 379                                 break;
 380                         q = isl_qpolynomial_mul(q, isl_qpolynomial_copy(base));
 381                         q = isl_qpolynomial_add(q, isl_qpolynomial_copy(r));
 382                 }
 383
 384                 if (isl_qpolynomial_is_zero(q))
 385                         sgn = isl_qpolynomial_sign(set, r);
 386                 else if (isl_qpolynomial_is_zero(r))
 387                         sgn = isl_qpolynomial_sign(set, q);
 388                 else {
 389                         int sgn_q, sgn_r;
 390                         sgn_r = isl_qpolynomial_sign(set, r);
 391                         sgn_q = isl_qpolynomial_sign(set, q);
 392                         if (sgn_r == sgn_q)
 393                                 sgn = sgn_r;
 394                 }
 395
 396                 isl_qpolynomial_free(min);
 397                 isl_qpolynomial_free(base);
 398                 isl_qpolynomial_free(q);
 399                 isl_qpolynomial_free(r);
 400         }
 401
 402         isl_int_clear(l);
 403
 404         isl_vec_free(v);
 405
 406         return sgn;
 407 }
 408
 409 /* Combine "fold1" and "fold2" into a single reduction, eliminating
 410  * those elements of one reduction that are already covered by the other
 411  * reduction on "set".
 412  *
 413  * If "fold1" or "fold2" is an empty reduction, then return
 414  * the other reduction.
 415  * If "fold1" or "fold2" is a NaN, then return this NaN.
 416  */
 417 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_fold_on_domain(
 418         __isl_keep isl_set *set,
 419         __isl_take isl_qpolynomial_fold *fold1,
 420         __isl_take isl_qpolynomial_fold *fold2)
 421 {
 422         int i, j;
 423         int n1;
 424         struct isl_qpolynomial_fold *res = NULL;
 425         int better;
 426
 427         if (!fold1 || !fold2)
 428                 goto error;
 429
 430         isl_assert(fold1->dim->ctx, fold1->type == fold2->type, goto error);
 431         isl_assert(fold1->dim->ctx, isl_space_is_equal(fold1->dim, fold2->dim),
 432                         goto error);
 433
 434         better = fold1->type == isl_fold_max ? -1 : 1;
 435
 436         if (isl_qpolynomial_fold_is_empty(fold1) ||
 437             isl_qpolynomial_fold_is_nan(fold2)) {
 438                 isl_qpolynomial_fold_free(fold1);
 439                 return fold2;
 440         }
 441
 442         if (isl_qpolynomial_fold_is_empty(fold2) ||
 443             isl_qpolynomial_fold_is_nan(fold1)) {
 444                 isl_qpolynomial_fold_free(fold2);
 445                 return fold1;
 446         }
 447
 448         res = qpolynomial_fold_alloc(fold1->type, isl_space_copy(fold1->dim),
 449                                         fold1->n + fold2->n);
 450         if (!res)
 451                 goto error;
 452
 453         for (i = 0; i < fold1->n; ++i) {
 454                 res->qp[res->n] = isl_qpolynomial_copy(fold1->qp[i]);
 455                 if (!res->qp[res->n])
 456                         goto error;
 457                 res->n++;
 458         }
 459         n1 = res->n;
 460
 461         for (i = 0; i < fold2->n; ++i) {
 462                 for (j = n1 - 1; j >= 0; --j) {
 463                         isl_qpolynomial *d;
 464                         int sgn, equal;
 465                         equal = isl_qpolynomial_plain_is_equal(res->qp[j],
 466                                                                 fold2->qp[i]);
 467                         if (equal < 0)
 468                                 goto error;
 469                         if (equal)
 470                                 break;
 471                         d = isl_qpolynomial_sub(
 472                                 isl_qpolynomial_copy(res->qp[j]),
 473                                 isl_qpolynomial_copy(fold2->qp[i]));
 474                         sgn = isl_qpolynomial_sign(set, d);
 475                         isl_qpolynomial_free(d);
 476                         if (sgn == 0)
 477                                 continue;
 478                         if (sgn != better)
 479                                 break;
 480                         isl_qpolynomial_free(res->qp[j]);
 481                         if (j != n1 - 1)
 482                                 res->qp[j] = res->qp[n1 - 1];
 483                         n1--;
 484                         if (n1 != res->n - 1)
 485                                 res->qp[n1] = res->qp[res->n - 1];
 486                         res->n--;
 487                 }
 488                 if (j >= 0)
 489                         continue;
 490                 res->qp[res->n] = isl_qpolynomial_copy(fold2->qp[i]);
 491                 if (!res->qp[res->n])
 492                         goto error;
 493                 res->n++;
 494         }
 495
 496         isl_qpolynomial_fold_free(fold1);
 497         isl_qpolynomial_fold_free(fold2);
 498
 499         return res;
 500 error:
 501         isl_qpolynomial_fold_free(res);
 502         isl_qpolynomial_fold_free(fold1);
 503         isl_qpolynomial_fold_free(fold2);
 504         return NULL;
 505 }
 506
 507 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_add_qpolynomial(
 508         __isl_take isl_qpolynomial_fold *fold, __isl_take isl_qpolynomial *qp)
 509 {
 510         int i;
 511
 512         if (!fold || !qp)
 513                 goto error;
 514
 515         if (isl_qpolynomial_is_zero(qp)) {
 516                 isl_qpolynomial_free(qp);
 517                 return fold;
 518         }
 519
 520         fold = isl_qpolynomial_fold_cow(fold);
 521         if (!fold)
 522                 goto error;
 523
 524         for (i = 0; i < fold->n; ++i) {
 525                 fold->qp[i] = isl_qpolynomial_add(fold->qp[i],
 526                                                 isl_qpolynomial_copy(qp));
 527                 if (!fold->qp[i])
 528                         goto error;
 529         }
 530
 531         isl_qpolynomial_free(qp);
 532         return fold;
 533 error:
 534         isl_qpolynomial_fold_free(fold);
 535         isl_qpolynomial_free(qp);
 536         return NULL;
 537 }
 538
 539 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_add_on_domain(
 540         __isl_keep isl_set *dom,
 541         __isl_take isl_qpolynomial_fold *fold1,
 542         __isl_take isl_qpolynomial_fold *fold2)
 543 {
 544         int i;
 545         isl_qpolynomial_fold *res = NULL;
 546
 547         if (!fold1 || !fold2)
 548                 goto error;
 549
 550         if (isl_qpolynomial_fold_is_empty(fold1)) {
 551                 isl_qpolynomial_fold_free(fold1);
 552                 return fold2;
 553         }
 554
 555         if (isl_qpolynomial_fold_is_empty(fold2)) {
 556                 isl_qpolynomial_fold_free(fold2);
 557                 return fold1;
 558         }
 559
 560         if (fold1->n == 1 && fold2->n != 1)
 561                 return isl_qpolynomial_fold_add_on_domain(dom, fold2, fold1);
 562
 563         if (fold2->n == 1) {
 564                 res = isl_qpolynomial_fold_add_qpolynomial(fold1,
 565                                         isl_qpolynomial_copy(fold2->qp[0]));
 566                 isl_qpolynomial_fold_free(fold2);
 567                 return res;
 568         }
 569
 570         res = isl_qpolynomial_fold_add_qpolynomial(
 571                                 isl_qpolynomial_fold_copy(fold1),
 572                                 isl_qpolynomial_copy(fold2->qp[0]));
 573
 574         for (i = 1; i < fold2->n; ++i) {
 575                 isl_qpolynomial_fold *res_i;
 576                 res_i = isl_qpolynomial_fold_add_qpolynomial(
 577                                         isl_qpolynomial_fold_copy(fold1),
 578                                         isl_qpolynomial_copy(fold2->qp[i]));
 579                 res = isl_qpolynomial_fold_fold_on_domain(dom, res, res_i);
 580         }
 581
 582         isl_qpolynomial_fold_free(fold1);
 583         isl_qpolynomial_fold_free(fold2);
 584         return res;
 585 error:
 586         isl_qpolynomial_fold_free(res);
 587         isl_qpolynomial_fold_free(fold1);
 588         isl_qpolynomial_fold_free(fold2);
 589         return NULL;
 590 }
 591
 592 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_substitute_equalities(
 593         __isl_take isl_qpolynomial_fold *fold, __isl_take isl_basic_set *eq)
 594 {
 595         int i;
 596
 597         if (!fold || !eq)
 598                 goto error;
 599
 600         fold = isl_qpolynomial_fold_cow(fold);
 601         if (!fold)
 602                 return NULL;
 603
 604         for (i = 0; i < fold->n; ++i) {
 605                 fold->qp[i] = isl_qpolynomial_substitute_equalities(fold->qp[i],
 606                                                         isl_basic_set_copy(eq));
 607                 if (!fold->qp[i])
 608                         goto error;
 609         }
 610
 611         isl_basic_set_free(eq);
 612         return fold;
 613 error:
 614         isl_basic_set_free(eq);
 615         isl_qpolynomial_fold_free(fold);
 616         return NULL;
 617 }
 618
 619 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_gist(
 620         __isl_take isl_qpolynomial_fold *fold, __isl_take isl_set *context)
 621 {
 622         int i;
 623
 624         if (!fold || !context)
 625                 goto error;
 626
 627         fold = isl_qpolynomial_fold_cow(fold);
 628         if (!fold)
 629                 return NULL;
 630
 631         for (i = 0; i < fold->n; ++i) {
 632                 fold->qp[i] = isl_qpolynomial_gist(fold->qp[i],
 633                                                         isl_set_copy(context));
 634                 if (!fold->qp[i])
 635                         goto error;
 636         }
 637
 638         isl_set_free(context);
 639         return fold;
 640 error:
 641         isl_set_free(context);
 642         isl_qpolynomial_fold_free(fold);
 643         return NULL;
 644 }
 645
 646 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_gist_params(
 647         __isl_take isl_qpolynomial_fold *fold, __isl_take isl_set *context)
 648 {
 649         isl_space *space = isl_qpolynomial_fold_get_domain_space(fold);
 650         isl_set *dom_context = isl_set_universe(space);
 651         dom_context = isl_set_intersect_params(dom_context, context);
 652         return isl_qpolynomial_fold_gist(fold, dom_context);
 653 }
 654
 655 #define HAS_TYPE
 656
 657 #undef PW
 658 #define PW isl_pw_qpolynomial_fold
 659 #undef EL
 660 #define EL isl_qpolynomial_fold
 661 #undef EL_IS_ZERO
 662 #define EL_IS_ZERO is_empty
 663 #undef ZERO
 664 #define ZERO zero
 665 #undef IS_ZERO
 666 #define IS_ZERO is_zero
 667 #undef FIELD
 668 #define FIELD fold
 669 #undef DEFAULT_IS_ZERO
 670 #define DEFAULT_IS_ZERO 1
 671
 672 #define NO_NEG
 673 #define NO_SUB
 674 #define NO_PULLBACK
 675
 676 #include <isl_pw_templ.c>
 677
 678 #undef UNION
 679 #define UNION isl_union_pw_qpolynomial_fold
 680 #undef PART
 681 #define PART isl_pw_qpolynomial_fold
 682 #undef PARTS
 683 #define PARTS pw_qpolynomial_fold
 684
 685 #define NO_SUB
 686
 687 #include <isl_union_single.c>
 688 #include <isl_union_eval.c>
 689
 690 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_empty(enum isl_fold type,
 691         __isl_take isl_space *dim)
 692 {
 693         return qpolynomial_fold_alloc(type, dim, 0);
 694 }
 695
 696 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_alloc(
 697         enum isl_fold type, __isl_take isl_qpolynomial *qp)
 698 {
 699         isl_qpolynomial_fold *fold;
 700
 701         if (!qp)
 702                 return NULL;
 703
 704         fold = qpolynomial_fold_alloc(type, isl_space_copy(qp->dim), 1);
 705         if (!fold)
 706                 goto error;
 707
 708         fold->qp[0] = qp;
 709         fold->n++;
 710
 711         return fold;
 712 error:
 713         isl_qpolynomial_fold_free(fold);
 714         isl_qpolynomial_free(qp);
 715         return NULL;
 716 }
 717
 718 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_copy(
 719         __isl_keep isl_qpolynomial_fold *fold)
 720 {
 721         if (!fold)
 722                 return NULL;
 723
 724         fold->ref++;
 725         return fold;
 726 }
 727
 728 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_dup(
 729         __isl_keep isl_qpolynomial_fold *fold)
 730 {
 731         int i;
 732         isl_qpolynomial_fold *dup;
 733
 734         if (!fold)
 735                 return NULL;
 736         dup = qpolynomial_fold_alloc(fold->type,
 737                                         isl_space_copy(fold->dim), fold->n);
 738         if (!dup)
 739                 return NULL;
 740
 741         dup->n = fold->n;
 742         for (i = 0; i < fold->n; ++i) {
 743                 dup->qp[i] = isl_qpolynomial_copy(fold->qp[i]);
 744                 if (!dup->qp[i])
 745                         goto error;
 746         }
 747
 748         return dup;
 749 error:
 750         isl_qpolynomial_fold_free(dup);
 751         return NULL;
 752 }
 753
 754 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_cow(
 755         __isl_take isl_qpolynomial_fold *fold)
 756 {
 757         if (!fold)
 758                 return NULL;
 759
 760         if (fold->ref == 1)
 761                 return fold;
 762         fold->ref--;
 763         return isl_qpolynomial_fold_dup(fold);
 764 }
 765
 766 void isl_qpolynomial_fold_free(__isl_take isl_qpolynomial_fold *fold)
 767 {
 768         int i;
 769
 770         if (!fold)
 771                 return;
 772         if (--fold->ref > 0)
 773                 return;
 774
 775         for (i = 0; i < fold->n; ++i)
 776                 isl_qpolynomial_free(fold->qp[i]);
 777         isl_space_free(fold->dim);
 778         free(fold);
 779 }
 780
 781 int isl_qpolynomial_fold_is_empty(__isl_keep isl_qpolynomial_fold *fold)
 782 {
 783         if (!fold)
 784                 return -1;
 785
 786         return fold->n == 0;
 787 }
 788
 789 /* Does "fold" represent max(NaN) or min(NaN)?
 790  */
 791 isl_bool isl_qpolynomial_fold_is_nan(__isl_keep isl_qpolynomial_fold *fold)
 792 {
 793         if (!fold)
 794                 return isl_bool_error;
 795         if (fold->n != 1)
 796                 return isl_bool_false;
 797         return isl_qpolynomial_is_nan(fold->qp[0]);
 798 }
 799
 800 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_fold(
 801         __isl_take isl_qpolynomial_fold *fold1,
 802         __isl_take isl_qpolynomial_fold *fold2)
 803 {
 804         int i;
 805         struct isl_qpolynomial_fold *res = NULL;
 806
 807         if (!fold1 || !fold2)
 808                 goto error;
 809
 810         isl_assert(fold1->dim->ctx, fold1->type == fold2->type, goto error);
 811         isl_assert(fold1->dim->ctx, isl_space_is_equal(fold1->dim, fold2->dim),
 812                         goto error);
 813
 814         if (isl_qpolynomial_fold_is_empty(fold1)) {
 815                 isl_qpolynomial_fold_free(fold1);
 816                 return fold2;
 817         }
 818
 819         if (isl_qpolynomial_fold_is_empty(fold2)) {
 820                 isl_qpolynomial_fold_free(fold2);
 821                 return fold1;
 822         }
 823
 824         res = qpolynomial_fold_alloc(fold1->type, isl_space_copy(fold1->dim),
 825                                         fold1->n + fold2->n);
 826         if (!res)
 827                 goto error;
 828
 829         for (i = 0; i < fold1->n; ++i) {
 830                 res->qp[res->n] = isl_qpolynomial_copy(fold1->qp[i]);
 831                 if (!res->qp[res->n])
 832                         goto error;
 833                 res->n++;
 834         }
 835
 836         for (i = 0; i < fold2->n; ++i) {
 837                 res->qp[res->n] = isl_qpolynomial_copy(fold2->qp[i]);
 838                 if (!res->qp[res->n])
 839                         goto error;
 840                 res->n++;
 841         }
 842
 843         isl_qpolynomial_fold_free(fold1);
 844         isl_qpolynomial_fold_free(fold2);
 845
 846         return res;
 847 error:
 848         isl_qpolynomial_fold_free(res);
 849         isl_qpolynomial_fold_free(fold1);
 850         isl_qpolynomial_fold_free(fold2);
 851         return NULL;
 852 }
 853
 854 __isl_give isl_pw_qpolynomial_fold *isl_pw_qpolynomial_fold_fold(
 855         __isl_take isl_pw_qpolynomial_fold *pw1,
 856         __isl_take isl_pw_qpolynomial_fold *pw2)
 857 {
 858         int i, j, n;
 859         struct isl_pw_qpolynomial_fold *res;
 860         isl_set *set;
 861
 862         if (!pw1 || !pw2)
 863                 goto error;
 864
 865         isl_assert(pw1->dim->ctx, isl_space_is_equal(pw1->dim, pw2->dim), goto error);
 866
 867         if (isl_pw_qpolynomial_fold_is_zero(pw1)) {
 868                 isl_pw_qpolynomial_fold_free(pw1);
 869                 return pw2;
 870         }
 871
 872         if (isl_pw_qpolynomial_fold_is_zero(pw2)) {
 873                 isl_pw_qpolynomial_fold_free(pw2);
 874                 return pw1;
 875         }
 876
 877         if (pw1->type != pw2->type)
 878                 isl_die(pw1->dim->ctx, isl_error_invalid,
 879                         "fold types don't match", goto error);
 880
 881         n = (pw1->n + 1) * (pw2->n + 1);
 882         res = isl_pw_qpolynomial_fold_alloc_size(isl_space_copy(pw1->dim),
 883                                                 pw1->type, n);
 884
 885         for (i = 0; i < pw1->n; ++i) {
 886                 set = isl_set_copy(pw1->p[i].set);
 887                 for (j = 0; j < pw2->n; ++j) {
 888                         struct isl_set *common;
 889                         isl_qpolynomial_fold *sum;
 890                         set = isl_set_subtract(set,
 891                                         isl_set_copy(pw2->p[j].set));
 892                         common = isl_set_intersect(isl_set_copy(pw1->p[i].set),
 893                                                 isl_set_copy(pw2->p[j].set));
 894                         if (isl_set_plain_is_empty(common)) {
 895                                 isl_set_free(common);
 896                                 continue;
 897                         }
 898
 899                         sum = isl_qpolynomial_fold_fold_on_domain(common,
 900                                isl_qpolynomial_fold_copy(pw1->p[i].fold),
 901                                isl_qpolynomial_fold_copy(pw2->p[j].fold));
 902
 903                         res = isl_pw_qpolynomial_fold_add_piece(res, common, sum);
 904                 }
 905                 res = isl_pw_qpolynomial_fold_add_piece(res, set,
 906                         isl_qpolynomial_fold_copy(pw1->p[i].fold));
 907         }
 908
 909         for (j = 0; j < pw2->n; ++j) {
 910                 set = isl_set_copy(pw2->p[j].set);
 911                 for (i = 0; i < pw1->n; ++i)
 912                         set = isl_set_subtract(set, isl_set_copy(pw1->p[i].set));
 913                 res = isl_pw_qpolynomial_fold_add_piece(res, set,
 914                                     isl_qpolynomial_fold_copy(pw2->p[j].fold));
 915         }
 916
 917         isl_pw_qpolynomial_fold_free(pw1);
 918         isl_pw_qpolynomial_fold_free(pw2);
 919
 920         return res;
 921 error:
 922         isl_pw_qpolynomial_fold_free(pw1);
 923         isl_pw_qpolynomial_fold_free(pw2);
 924         return NULL;
 925 }
 926
 927 __isl_give isl_union_pw_qpolynomial_fold *isl_union_pw_qpolynomial_fold_fold_pw_qpolynomial_fold(
 928         __isl_take isl_union_pw_qpolynomial_fold *u,
 929         __isl_take isl_pw_qpolynomial_fold *part)
 930 {
 931         struct isl_hash_table_entry *entry;
 932
 933         u = isl_union_pw_qpolynomial_fold_cow(u);
 934
 935         if (!part || !u)
 936                 goto error;
 937
 938         isl_assert(u->space->ctx,
 939             isl_space_match(part->dim, isl_dim_param, u->space, isl_dim_param),
 940             goto error);
 941
 942         entry = isl_union_pw_qpolynomial_fold_find_part_entry(u, part->dim, 1);
 943         if (!entry)
 944                 goto error;
 945
 946         if (!entry->data)
 947                 entry->data = part;
 948         else {
 949                 entry->data = isl_pw_qpolynomial_fold_fold(entry->data,
 950                                             isl_pw_qpolynomial_fold_copy(part));
 951                 if (!entry->data)
 952                         goto error;
 953                 isl_pw_qpolynomial_fold_free(part);
 954         }
 955
 956         return u;
 957 error:
 958         isl_pw_qpolynomial_fold_free(part);
 959         isl_union_pw_qpolynomial_fold_free(u);
 960         return NULL;
 961 }
 962
 963 static isl_stat fold_part(__isl_take isl_pw_qpolynomial_fold *part, void *user)
 964 {
 965         isl_union_pw_qpolynomial_fold **u;
 966         u = (isl_union_pw_qpolynomial_fold **)user;
 967
 968         *u = isl_union_pw_qpolynomial_fold_fold_pw_qpolynomial_fold(*u, part);
 969
 970         return isl_stat_ok;
 971 }
 972
 973 __isl_give isl_union_pw_qpolynomial_fold *isl_union_pw_qpolynomial_fold_fold(
 974         __isl_take isl_union_pw_qpolynomial_fold *u1,
 975         __isl_take isl_union_pw_qpolynomial_fold *u2)
 976 {
 977         u1 = isl_union_pw_qpolynomial_fold_cow(u1);
 978
 979         if (!u1 || !u2)
 980                 goto error;
 981
 982         if (isl_union_pw_qpolynomial_fold_foreach_pw_qpolynomial_fold(u2,
 983                                                         &fold_part, &u1) < 0)
 984                 goto error;
 985
 986         isl_union_pw_qpolynomial_fold_free(u2);
 987
 988         return u1;
 989 error:
 990         isl_union_pw_qpolynomial_fold_free(u1);
 991         isl_union_pw_qpolynomial_fold_free(u2);
 992         return NULL;
 993 }
 994
 995 __isl_give isl_pw_qpolynomial_fold *isl_pw_qpolynomial_fold_from_pw_qpolynomial(
 996         enum isl_fold type, __isl_take isl_pw_qpolynomial *pwqp)
 997 {
 998         int i;
 999         isl_pw_qpolynomial_fold *pwf;
1000
1001         if (!pwqp)
1002                 return NULL;
1003
1004         pwf = isl_pw_qpolynomial_fold_alloc_size(isl_space_copy(pwqp->dim),
1005                                                 type, pwqp->n);
1006
1007         for (i = 0; i < pwqp->n; ++i)
1008                 pwf = isl_pw_qpolynomial_fold_add_piece(pwf,
1009                         isl_set_copy(pwqp->p[i].set),
1010                         isl_qpolynomial_fold_alloc(type,
1011                                 isl_qpolynomial_copy(pwqp->p[i].qp)));
1012
1013         isl_pw_qpolynomial_free(pwqp);
1014
1015         return pwf;
1016 }
1017
1018 __isl_give isl_pw_qpolynomial_fold *isl_pw_qpolynomial_fold_add(
1019         __isl_take isl_pw_qpolynomial_fold *pwf1,
1020         __isl_take isl_pw_qpolynomial_fold *pwf2)
1021 {
1022         return isl_pw_qpolynomial_fold_union_add_(pwf1, pwf2);
1023 }
1024
1025 int isl_qpolynomial_fold_plain_is_equal(__isl_keep isl_qpolynomial_fold *fold1,
1026         __isl_keep isl_qpolynomial_fold *fold2)
1027 {
1028         int i;
1029
1030         if (!fold1 || !fold2)
1031                 return -1;
1032
1033         if (fold1->n != fold2->n)
1034                 return 0;
1035
1036         /* We probably want to sort the qps first... */
1037         for (i = 0; i < fold1->n; ++i) {
1038                 int eq = isl_qpolynomial_plain_is_equal(fold1->qp[i], fold2->qp[i]);
1039                 if (eq < 0 || !eq)
1040                         return eq;
1041         }
1042
1043         return 1;
1044 }
1045
1046 __isl_give isl_val *isl_qpolynomial_fold_eval(
1047         __isl_take isl_qpolynomial_fold *fold, __isl_take isl_point *pnt)
1048 {
1049         isl_ctx *ctx;
1050         isl_val *v;
1051
1052         if (!fold || !pnt)
1053                 goto error;
1054         ctx = isl_point_get_ctx(pnt);
1055         isl_assert(pnt->dim->ctx, isl_space_is_equal(pnt->dim, fold->dim), goto error);
1056         isl_assert(pnt->dim->ctx,
1057                 fold->type == isl_fold_max || fold->type == isl_fold_min,
1058                 goto error);
1059
1060         if (fold->n == 0)
1061                 v = isl_val_zero(ctx);
1062         else {
1063                 int i;
1064                 v = isl_qpolynomial_eval(isl_qpolynomial_copy(fold->qp[0]),
1065                                                 isl_point_copy(pnt));
1066                 for (i = 1; i < fold->n; ++i) {
1067                         isl_val *v_i;
1068                         v_i = isl_qpolynomial_eval(
1069                                             isl_qpolynomial_copy(fold->qp[i]),
1070                                             isl_point_copy(pnt));
1071                         if (fold->type == isl_fold_max)
1072                                 v = isl_val_max(v, v_i);
1073                         else
1074                                 v = isl_val_min(v, v_i);
1075                 }
1076         }
1077         isl_qpolynomial_fold_free(fold);
1078         isl_point_free(pnt);
1079
1080         return v;
1081 error:
1082         isl_qpolynomial_fold_free(fold);
1083         isl_point_free(pnt);
1084         return NULL;
1085 }
1086
1087 size_t isl_pw_qpolynomial_fold_size(__isl_keep isl_pw_qpolynomial_fold *pwf)
1088 {
1089         int i;
1090         size_t n = 0;
1091
1092         for (i = 0; i < pwf->n; ++i)
1093                 n += pwf->p[i].fold->n;
1094
1095         return n;
1096 }
1097
1098 __isl_give isl_val *isl_qpolynomial_fold_opt_on_domain(
1099         __isl_take isl_qpolynomial_fold *fold, __isl_take isl_set *set, int max)
1100 {
1101         int i;
1102         isl_val *opt;
1103
1104         if (!set || !fold)
1105                 goto error;
1106
1107         if (fold->n == 0) {
1108                 opt = isl_val_zero(isl_set_get_ctx(set));
1109                 isl_set_free(set);
1110                 isl_qpolynomial_fold_free(fold);
1111                 return opt;
1112         }
1113
1114         opt = isl_qpolynomial_opt_on_domain(isl_qpolynomial_copy(fold->qp[0]),
1115                                                 isl_set_copy(set), max);
1116         for (i = 1; i < fold->n; ++i) {
1117                 isl_val *opt_i;
1118                 opt_i = isl_qpolynomial_opt_on_domain(
1119                                 isl_qpolynomial_copy(fold->qp[i]),
1120                                 isl_set_copy(set), max);
1121                 if (max)
1122                         opt = isl_val_max(opt, opt_i);
1123                 else
1124                         opt = isl_val_min(opt, opt_i);
1125         }
1126
1127         isl_set_free(set);
1128         isl_qpolynomial_fold_free(fold);
1129
1130         return opt;
1131 error:
1132         isl_set_free(set);
1133         isl_qpolynomial_fold_free(fold);
1134         return NULL;
1135 }
1136
1137 /* Check whether for each quasi-polynomial in "fold2" there is
1138  * a quasi-polynomial in "fold1" that dominates it on "set".
1139  */
1140 static int qpolynomial_fold_covers_on_domain(__isl_keep isl_set *set,
1141         __isl_keep isl_qpolynomial_fold *fold1,
1142         __isl_keep isl_qpolynomial_fold *fold2)
1143 {
1144         int i, j;
1145         int covers;
1146
1147         if (!set || !fold1 || !fold2)
1148                 return -1;
1149
1150         covers = fold1->type == isl_fold_max ? 1 : -1;
1151
1152         for (i = 0; i < fold2->n; ++i) {
1153                 for (j = 0; j < fold1->n; ++j) {
1154                         isl_qpolynomial *d;
1155                         int sgn;
1156
1157                         d = isl_qpolynomial_sub(
1158                                 isl_qpolynomial_copy(fold1->qp[j]),
1159                                 isl_qpolynomial_copy(fold2->qp[i]));
1160                         sgn = isl_qpolynomial_sign(set, d);
1161                         isl_qpolynomial_free(d);
1162                         if (sgn == covers)
1163                                 break;
1164                 }
1165                 if (j >= fold1->n)
1166                         return 0;
1167         }
1168
1169         return 1;
1170 }
1171
1172 /* Check whether "pwf1" dominated "pwf2", i.e., the domain of "pwf1" contains
1173  * that of "pwf2" and on each cell, the corresponding fold from pwf1 dominates
1174  * that of pwf2.
1175  */
1176 int isl_pw_qpolynomial_fold_covers(__isl_keep isl_pw_qpolynomial_fold *pwf1,
1177         __isl_keep isl_pw_qpolynomial_fold *pwf2)
1178 {
1179         int i, j;
1180         isl_set *dom1, *dom2;
1181         int is_subset;
1182
1183         if (!pwf1 || !pwf2)
1184                 return -1;
1185
1186         if (pwf2->n == 0)
1187                 return 1;
1188         if (pwf1->n == 0)
1189                 return 0;
1190
1191         dom1 = isl_pw_qpolynomial_fold_domain(isl_pw_qpolynomial_fold_copy(pwf1));
1192         dom2 = isl_pw_qpolynomial_fold_domain(isl_pw_qpolynomial_fold_copy(pwf2));
1193         is_subset = isl_set_is_subset(dom2, dom1);
1194         isl_set_free(dom1);
1195         isl_set_free(dom2);
1196
1197         if (is_subset < 0 || !is_subset)
1198                 return is_subset;
1199
1200         for (i = 0; i < pwf2->n; ++i) {
1201                 for (j = 0; j < pwf1->n; ++j) {
1202                         int is_empty;
1203                         isl_set *common;
1204                         int covers;
1205
1206                         common = isl_set_intersect(isl_set_copy(pwf1->p[j].set),
1207                                                    isl_set_copy(pwf2->p[i].set));
1208                         is_empty = isl_set_is_empty(common);
1209                         if (is_empty < 0 || is_empty) {
1210                                 isl_set_free(common);
1211                                 if (is_empty < 0)
1212                                         return -1;
1213                                 continue;
1214                         }
1215                         covers = qpolynomial_fold_covers_on_domain(common,
1216                                         pwf1->p[j].fold, pwf2->p[i].fold);
1217                         isl_set_free(common);
1218                         if (covers < 0 || !covers)
1219                                 return covers;
1220                 }
1221         }
1222
1223         return 1;
1224 }
1225
1226 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_morph_domain(
1227         __isl_take isl_qpolynomial_fold *fold, __isl_take isl_morph *morph)
1228 {
1229         int i;
1230         isl_ctx *ctx;
1231
1232         if (!fold || !morph)
1233                 goto error;
1234
1235         ctx = fold->dim->ctx;
1236         isl_assert(ctx, isl_space_is_equal(fold->dim, morph->dom->dim), goto error);
1237
1238         fold = isl_qpolynomial_fold_cow(fold);
1239         if (!fold)
1240                 goto error;
1241
1242         isl_space_free(fold->dim);
1243         fold->dim = isl_space_copy(morph->ran->dim);
1244         if (!fold->dim)
1245                 goto error;
1246
1247         for (i = 0; i < fold->n; ++i) {
1248                 fold->qp[i] = isl_qpolynomial_morph_domain(fold->qp[i],
1249                                                 isl_morph_copy(morph));
1250                 if (!fold->qp[i])
1251                         goto error;
1252         }
1253
1254         isl_morph_free(morph);
1255
1256         return fold;
1257 error:
1258         isl_qpolynomial_fold_free(fold);
1259         isl_morph_free(morph);
1260         return NULL;
1261 }
1262
1263 enum isl_fold isl_qpolynomial_fold_get_type(__isl_keep isl_qpolynomial_fold *fold)
1264 {
1265         if (!fold)
1266                 return isl_fold_list;
1267         return fold->type;
1268 }
1269
1270 enum isl_fold isl_union_pw_qpolynomial_fold_get_type(
1271         __isl_keep isl_union_pw_qpolynomial_fold *upwf)
1272 {
1273         if (!upwf)
1274                 return isl_fold_list;
1275         return upwf->type;
1276 }
1277
1278 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_lift(
1279         __isl_take isl_qpolynomial_fold *fold, __isl_take isl_space *dim)
1280 {
1281         int i;
1282
1283         if (!fold || !dim)
1284                 goto error;
1285
1286         if (isl_space_is_equal(fold->dim, dim)) {
1287                 isl_space_free(dim);
1288                 return fold;
1289         }
1290
1291         fold = isl_qpolynomial_fold_cow(fold);
1292         if (!fold)
1293                 goto error;
1294
1295         isl_space_free(fold->dim);
1296         fold->dim = isl_space_copy(dim);
1297         if (!fold->dim)
1298                 goto error;
1299
1300         for (i = 0; i < fold->n; ++i) {
1301                 fold->qp[i] = isl_qpolynomial_lift(fold->qp[i],
1302                                                 isl_space_copy(dim));
1303                 if (!fold->qp[i])
1304                         goto error;
1305         }
1306
1307         isl_space_free(dim);
1308
1309         return fold;
1310 error:
1311         isl_qpolynomial_fold_free(fold);
1312         isl_space_free(dim);
1313         return NULL;
1314 }
1315
1316 isl_stat isl_qpolynomial_fold_foreach_qpolynomial(
1317         __isl_keep isl_qpolynomial_fold *fold,
1318         isl_stat (*fn)(__isl_take isl_qpolynomial *qp, void *user), void *user)
1319 {
1320         int i;
1321
1322         if (!fold)
1323                 return isl_stat_error;
1324
1325         for (i = 0; i < fold->n; ++i)
1326                 if (fn(isl_qpolynomial_copy(fold->qp[i]), user) < 0)
1327                         return isl_stat_error;
1328
1329         return isl_stat_ok;
1330 }
1331
1332 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_move_dims(
1333         __isl_take isl_qpolynomial_fold *fold,
1334         enum isl_dim_type dst_type, unsigned dst_pos,
1335         enum isl_dim_type src_type, unsigned src_pos, unsigned n)
1336 {
1337         int i;
1338
1339         if (n == 0)
1340                 return fold;
1341
1342         fold = isl_qpolynomial_fold_cow(fold);
1343         if (!fold)
1344                 return NULL;
1345
1346         fold->dim = isl_space_move_dims(fold->dim, dst_type, dst_pos,
1347                                                 src_type, src_pos, n);
1348         if (!fold->dim)
1349                 goto error;
1350
1351         for (i = 0; i < fold->n; ++i) {
1352                 fold->qp[i] = isl_qpolynomial_move_dims(fold->qp[i],
1353                                 dst_type, dst_pos, src_type, src_pos, n);
1354                 if (!fold->qp[i])
1355                         goto error;
1356         }
1357
1358         return fold;
1359 error:
1360         isl_qpolynomial_fold_free(fold);
1361         return NULL;
1362 }
1363
1364 /* For each 0 <= i < "n", replace variable "first" + i of type "type"
1365  * in fold->qp[k] by subs[i].
1366  */
1367 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_substitute(
1368         __isl_take isl_qpolynomial_fold *fold,
1369         enum isl_dim_type type, unsigned first, unsigned n,
1370         __isl_keep isl_qpolynomial **subs)
1371 {
1372         int i;
1373
1374         if (n == 0)
1375                 return fold;
1376
1377         fold = isl_qpolynomial_fold_cow(fold);
1378         if (!fold)
1379                 return NULL;
1380
1381         for (i = 0; i < fold->n; ++i) {
1382                 fold->qp[i] = isl_qpolynomial_substitute(fold->qp[i],
1383                                 type, first, n, subs);
1384                 if (!fold->qp[i])
1385                         goto error;
1386         }
1387
1388         return fold;
1389 error:
1390         isl_qpolynomial_fold_free(fold);
1391         return NULL;
1392 }
1393
1394 static isl_stat add_pwqp(__isl_take isl_pw_qpolynomial *pwqp, void *user)
1395 {
1396         isl_ctx *ctx;
1397         isl_pw_qpolynomial_fold *pwf;
1398         isl_union_pw_qpolynomial_fold **upwf;
1399         struct isl_hash_table_entry *entry;
1400
1401         upwf = (isl_union_pw_qpolynomial_fold **)user;
1402
1403         ctx = pwqp->dim->ctx;
1404         entry = isl_union_pw_qpolynomial_fold_find_part_entry(*upwf,
1405                          pwqp->dim, 1);
1406         if (!entry)
1407                 goto error;
1408
1409         pwf = isl_pw_qpolynomial_fold_from_pw_qpolynomial((*upwf)->type, pwqp);
1410         if (!entry->data)
1411                 entry->data = pwf;
1412         else {
1413                 entry->data = isl_pw_qpolynomial_fold_add(entry->data, pwf);
1414                 if (!entry->data)
1415                         return isl_stat_error;
1416                 if (isl_pw_qpolynomial_fold_is_zero(entry->data))
1417                         *upwf = isl_union_pw_qpolynomial_fold_remove_part_entry(
1418                                                                 *upwf, entry);
1419         }
1420
1421         return isl_stat_ok;
1422 error:
1423         isl_pw_qpolynomial_free(pwqp);
1424         return isl_stat_error;
1425 }
1426
1427 __isl_give isl_union_pw_qpolynomial_fold *isl_union_pw_qpolynomial_fold_add_union_pw_qpolynomial(
1428         __isl_take isl_union_pw_qpolynomial_fold *upwf,
1429         __isl_take isl_union_pw_qpolynomial *upwqp)
1430 {
1431         upwf = isl_union_pw_qpolynomial_fold_align_params(upwf,
1432                                 isl_union_pw_qpolynomial_get_space(upwqp));
1433         upwqp = isl_union_pw_qpolynomial_align_params(upwqp,
1434                                 isl_union_pw_qpolynomial_fold_get_space(upwf));
1435
1436         upwf = isl_union_pw_qpolynomial_fold_cow(upwf);
1437         if (!upwf || !upwqp)
1438                 goto error;
1439
1440         if (isl_union_pw_qpolynomial_foreach_pw_qpolynomial(upwqp, &add_pwqp,
1441                                                          &upwf) < 0)
1442                 goto error;
1443
1444         isl_union_pw_qpolynomial_free(upwqp);
1445
1446         return upwf;
1447 error:
1448         isl_union_pw_qpolynomial_fold_free(upwf);
1449         isl_union_pw_qpolynomial_free(upwqp);
1450         return NULL;
1451 }
1452
1453 static int join_compatible(__isl_keep isl_space *dim1, __isl_keep isl_space *dim2)
1454 {
1455         int m;
1456         m = isl_space_match(dim1, isl_dim_param, dim2, isl_dim_param);
1457         if (m < 0 || !m)
1458                 return m;
1459         return isl_space_tuple_is_equal(dim1, isl_dim_out, dim2, isl_dim_in);
1460 }
1461
1462 /* Compute the intersection of the range of the map and the domain
1463  * of the piecewise quasipolynomial reduction and then compute a bound
1464  * on the associated quasipolynomial reduction over all elements
1465  * in this intersection.
1466  *
1467  * We first introduce some unconstrained dimensions in the
1468  * piecewise quasipolynomial, intersect the resulting domain
1469  * with the wrapped map and the compute the sum.
1470  */
1471 __isl_give isl_pw_qpolynomial_fold *isl_map_apply_pw_qpolynomial_fold(
1472         __isl_take isl_map *map, __isl_take isl_pw_qpolynomial_fold *pwf,
1473         int *tight)
1474 {
1475         isl_ctx *ctx;
1476         isl_set *dom;
1477         isl_space *map_dim;
1478         isl_space *pwf_dim;
1479         unsigned n_in;
1480         int ok;
1481
1482         ctx = isl_map_get_ctx(map);
1483         if (!ctx)
1484                 goto error;
1485
1486         map_dim = isl_map_get_space(map);
1487         pwf_dim = isl_pw_qpolynomial_fold_get_space(pwf);
1488         ok = join_compatible(map_dim, pwf_dim);
1489         isl_space_free(map_dim);
1490         isl_space_free(pwf_dim);
1491         if (!ok)
1492                 isl_die(ctx, isl_error_invalid, "incompatible dimensions",
1493                         goto error);
1494
1495         n_in = isl_map_dim(map, isl_dim_in);
1496         pwf = isl_pw_qpolynomial_fold_insert_dims(pwf, isl_dim_in, 0, n_in);
1497
1498         dom = isl_map_wrap(map);
1499         pwf = isl_pw_qpolynomial_fold_reset_domain_space(pwf,
1500                                                 isl_set_get_space(dom));
1501
1502         pwf = isl_pw_qpolynomial_fold_intersect_domain(pwf, dom);
1503         pwf = isl_pw_qpolynomial_fold_bound(pwf, tight);
1504
1505         return pwf;
1506 error:
1507         isl_map_free(map);
1508         isl_pw_qpolynomial_fold_free(pwf);
1509         return NULL;
1510 }
1511
1512 __isl_give isl_pw_qpolynomial_fold *isl_set_apply_pw_qpolynomial_fold(
1513         __isl_take isl_set *set, __isl_take isl_pw_qpolynomial_fold *pwf,
1514         int *tight)
1515 {
1516         return isl_map_apply_pw_qpolynomial_fold(set, pwf, tight);
1517 }
1518
1519 struct isl_apply_fold_data {
1520         isl_union_pw_qpolynomial_fold *upwf;
1521         isl_union_pw_qpolynomial_fold *res;
1522         isl_map *map;
1523         int tight;
1524 };
1525
1526 static isl_stat pw_qpolynomial_fold_apply(
1527         __isl_take isl_pw_qpolynomial_fold *pwf, void *user)
1528 {
1529         isl_space *map_dim;
1530         isl_space *pwf_dim;
1531         struct isl_apply_fold_data *data = user;
1532         int ok;
1533
1534         map_dim = isl_map_get_space(data->map);
1535         pwf_dim = isl_pw_qpolynomial_fold_get_space(pwf);
1536         ok = join_compatible(map_dim, pwf_dim);
1537         isl_space_free(map_dim);
1538         isl_space_free(pwf_dim);
1539
1540         if (ok) {
1541                 pwf = isl_map_apply_pw_qpolynomial_fold(isl_map_copy(data->map),
1542                                     pwf, data->tight ? &data->tight : NULL);
1543                 data->res = isl_union_pw_qpolynomial_fold_fold_pw_qpolynomial_fold(
1544                                                         data->res, pwf);
1545         } else
1546                 isl_pw_qpolynomial_fold_free(pwf);
1547
1548         return isl_stat_ok;
1549 }
1550
1551 static isl_stat map_apply(__isl_take isl_map *map, void *user)
1552 {
1553         struct isl_apply_fold_data *data = user;
1554         isl_stat r;
1555
1556         data->map = map;
1557         r = isl_union_pw_qpolynomial_fold_foreach_pw_qpolynomial_fold(
1558                                 data->upwf, &pw_qpolynomial_fold_apply, data);
1559
1560         isl_map_free(map);
1561         return r;
1562 }
1563
1564 __isl_give isl_union_pw_qpolynomial_fold *isl_union_map_apply_union_pw_qpolynomial_fold(
1565         __isl_take isl_union_map *umap,
1566         __isl_take isl_union_pw_qpolynomial_fold *upwf, int *tight)
1567 {
1568         isl_space *dim;
1569         enum isl_fold type;
1570         struct isl_apply_fold_data data;
1571
1572         upwf = isl_union_pw_qpolynomial_fold_align_params(upwf,
1573                                 isl_union_map_get_space(umap));
1574         umap = isl_union_map_align_params(umap,
1575                                 isl_union_pw_qpolynomial_fold_get_space(upwf));
1576
1577         data.upwf = upwf;
1578         data.tight = tight ? 1 : 0;
1579         dim = isl_union_pw_qpolynomial_fold_get_space(upwf);
1580         type = isl_union_pw_qpolynomial_fold_get_type(upwf);
1581         data.res = isl_union_pw_qpolynomial_fold_zero(dim, type);
1582         if (isl_union_map_foreach_map(umap, &map_apply, &data) < 0)
1583                 goto error;
1584
1585         isl_union_map_free(umap);
1586         isl_union_pw_qpolynomial_fold_free(upwf);
1587
1588         if (tight)
1589                 *tight = data.tight;
1590
1591         return data.res;
1592 error:
1593         isl_union_map_free(umap);
1594         isl_union_pw_qpolynomial_fold_free(upwf);
1595         isl_union_pw_qpolynomial_fold_free(data.res);
1596         return NULL;
1597 }
1598
1599 __isl_give isl_union_pw_qpolynomial_fold *isl_union_set_apply_union_pw_qpolynomial_fold(
1600         __isl_take isl_union_set *uset,
1601         __isl_take isl_union_pw_qpolynomial_fold *upwf, int *tight)
1602 {
1603         return isl_union_map_apply_union_pw_qpolynomial_fold(uset, upwf, tight);
1604 }
1605
1606 /* Reorder the dimension of "fold" according to the given reordering.
1607  */
1608 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_realign_domain(
1609         __isl_take isl_qpolynomial_fold *fold, __isl_take isl_reordering *r)
1610 {
1611         int i;
1612
1613         fold = isl_qpolynomial_fold_cow(fold);
1614         if (!fold || !r)
1615                 goto error;
1616
1617         for (i = 0; i < fold->n; ++i) {
1618                 fold->qp[i] = isl_qpolynomial_realign_domain(fold->qp[i],
1619                                                     isl_reordering_copy(r));
1620                 if (!fold->qp[i])
1621                         goto error;
1622         }
1623
1624         fold = isl_qpolynomial_fold_reset_domain_space(fold,
1625                                                     isl_space_copy(r->dim));
1626
1627         isl_reordering_free(r);
1628
1629         return fold;
1630 error:
1631         isl_qpolynomial_fold_free(fold);
1632         isl_reordering_free(r);
1633         return NULL;
1634 }
1635
1636 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_mul_isl_int(
1637         __isl_take isl_qpolynomial_fold *fold, isl_int v)
1638 {
1639         int i;
1640
1641         if (isl_int_is_one(v))
1642                 return fold;
1643         if (fold && isl_int_is_zero(v)) {
1644                 isl_qpolynomial_fold *zero;
1645                 isl_space *dim = isl_space_copy(fold->dim);
1646                 zero = isl_qpolynomial_fold_empty(fold->type, dim);
1647                 isl_qpolynomial_fold_free(fold);
1648                 return zero;
1649         }
1650
1651         fold = isl_qpolynomial_fold_cow(fold);
1652         if (!fold)
1653                 return NULL;
1654
1655         if (isl_int_is_neg(v))
1656                 fold->type = isl_fold_type_negate(fold->type);
1657         for (i = 0; i < fold->n; ++i) {
1658                 fold->qp[i] = isl_qpolynomial_mul_isl_int(fold->qp[i], v);
1659                 if (!fold->qp[i])
1660                         goto error;
1661         }
1662
1663         return fold;
1664 error:
1665         isl_qpolynomial_fold_free(fold);
1666         return NULL;
1667 }
1668
1669 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_scale(
1670         __isl_take isl_qpolynomial_fold *fold, isl_int v)
1671 {
1672         return isl_qpolynomial_fold_mul_isl_int(fold, v);
1673 }
1674
1675 /* Multiply "fold" by "v".
1676  */
1677 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_scale_val(
1678         __isl_take isl_qpolynomial_fold *fold, __isl_take isl_val *v)
1679 {
1680         int i;
1681
1682         if (!fold || !v)
1683                 goto error;
1684
1685         if (isl_val_is_one(v)) {
1686                 isl_val_free(v);
1687                 return fold;
1688         }
1689         if (isl_val_is_zero(v)) {
1690                 isl_qpolynomial_fold *zero;
1691                 isl_space *space = isl_qpolynomial_fold_get_domain_space(fold);
1692                 zero = isl_qpolynomial_fold_empty(fold->type, space);
1693                 isl_qpolynomial_fold_free(fold);
1694                 isl_val_free(v);
1695                 return zero;
1696         }
1697         if (!isl_val_is_rat(v))
1698                 isl_die(isl_qpolynomial_fold_get_ctx(fold), isl_error_invalid,
1699                         "expecting rational factor", goto error);
1700
1701         fold = isl_qpolynomial_fold_cow(fold);
1702         if (!fold)
1703                 goto error;
1704
1705         if (isl_val_is_neg(v))
1706                 fold->type = isl_fold_type_negate(fold->type);
1707         for (i = 0; i < fold->n; ++i) {
1708                 fold->qp[i] = isl_qpolynomial_scale_val(fold->qp[i],
1709                                                         isl_val_copy(v));
1710                 if (!fold->qp[i])
1711                         goto error;
1712         }
1713
1714         isl_val_free(v);
1715         return fold;
1716 error:
1717         isl_val_free(v);
1718         isl_qpolynomial_fold_free(fold);
1719         return NULL;
1720 }
1721
1722 /* Divide "fold" by "v".
1723  */
1724 __isl_give isl_qpolynomial_fold *isl_qpolynomial_fold_scale_down_val(
1725         __isl_take isl_qpolynomial_fold *fold, __isl_take isl_val *v)
1726 {
1727         if (!fold || !v)
1728                 goto error;
1729
1730         if (isl_val_is_one(v)) {
1731                 isl_val_free(v);
1732                 return fold;
1733         }
1734         if (!isl_val_is_rat(v))
1735                 isl_die(isl_qpolynomial_fold_get_ctx(fold), isl_error_invalid,
1736                         "expecting rational factor", goto error);
1737         if (isl_val_is_zero(v))
1738                 isl_die(isl_val_get_ctx(v), isl_error_invalid,
1739                         "cannot scale down by zero", goto error);
1740
1741         return isl_qpolynomial_fold_scale_val(fold, isl_val_inv(v));
1742 error:
1743         isl_val_free(v);
1744         isl_qpolynomial_fold_free(fold);
1745         return NULL;
1746 }