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