isl_morph.c

   1 /*
   2  * Copyright 2010-2011 INRIA Saclay
   3  * Copyright 2014      Ecole Normale Superieure
   4  *
   5  * Use of this software is governed by the MIT license
   6  *
   7  * Written by Sven Verdoolaege, INRIA Saclay - Ile-de-France,
   8  * Parc Club Orsay Universite, ZAC des vignes, 4 rue Jacques Monod,
   9  * 91893 Orsay, France
  10  * and Ecole Normale Superieure, 45 rue d'Ulm, 75230 Paris, France
  11  */
  12
  13 #include <isl_map_private.h>
  14 #include <isl_aff_private.h>
  15 #include <isl_morph.h>
  16 #include <isl_seq.h>
  17 #include <isl_mat_private.h>
  18 #include <isl_space_private.h>
  19 #include <isl_equalities.h>
  20
  21 isl_ctx *isl_morph_get_ctx(__isl_keep isl_morph *morph)
  22 {
  23         if (!morph)
  24                 return NULL;
  25         return isl_basic_set_get_ctx(morph->dom);
  26 }
  27
  28 __isl_give isl_morph *isl_morph_alloc(
  29         __isl_take isl_basic_set *dom, __isl_take isl_basic_set *ran,
  30         __isl_take isl_mat *map, __isl_take isl_mat *inv)
  31 {
  32         isl_morph *morph;
  33
  34         if (!dom || !ran || !map || !inv)
  35                 goto error;
  36
  37         morph = isl_alloc_type(dom->ctx, struct isl_morph);
  38         if (!morph)
  39                 goto error;
  40
  41         morph->ref = 1;
  42         morph->dom = dom;
  43         morph->ran = ran;
  44         morph->map = map;
  45         morph->inv = inv;
  46
  47         return morph;
  48 error:
  49         isl_basic_set_free(dom);
  50         isl_basic_set_free(ran);
  51         isl_mat_free(map);
  52         isl_mat_free(inv);
  53         return NULL;
  54 }
  55
  56 __isl_give isl_morph *isl_morph_copy(__isl_keep isl_morph *morph)
  57 {
  58         if (!morph)
  59                 return NULL;
  60
  61         morph->ref++;
  62         return morph;
  63 }
  64
  65 __isl_give isl_morph *isl_morph_dup(__isl_keep isl_morph *morph)
  66 {
  67         if (!morph)
  68                 return NULL;
  69
  70         return isl_morph_alloc(isl_basic_set_copy(morph->dom),
  71                 isl_basic_set_copy(morph->ran),
  72                 isl_mat_copy(morph->map), isl_mat_copy(morph->inv));
  73 }
  74
  75 __isl_give isl_morph *isl_morph_cow(__isl_take isl_morph *morph)
  76 {
  77         if (!morph)
  78                 return NULL;
  79
  80         if (morph->ref == 1)
  81                 return morph;
  82         morph->ref--;
  83         return isl_morph_dup(morph);
  84 }
  85
  86 __isl_null isl_morph *isl_morph_free(__isl_take isl_morph *morph)
  87 {
  88         if (!morph)
  89                 return NULL;
  90
  91         if (--morph->ref > 0)
  92                 return NULL;
  93
  94         isl_basic_set_free(morph->dom);
  95         isl_basic_set_free(morph->ran);
  96         isl_mat_free(morph->map);
  97         isl_mat_free(morph->inv);
  98         free(morph);
  99
 100         return NULL;
 101 }
 102
 103 /* Is "morph" an identity on the parameters?
 104  */
 105 static int identity_on_parameters(__isl_keep isl_morph *morph)
 106 {
 107         int is_identity;
 108         unsigned nparam;
 109         isl_mat *sub;
 110
 111         nparam = isl_morph_dom_dim(morph, isl_dim_param);
 112         if (nparam != isl_morph_ran_dim(morph, isl_dim_param))
 113                 return 0;
 114         if (nparam == 0)
 115                 return 1;
 116         sub = isl_mat_sub_alloc(morph->map, 0, 1 + nparam, 0, 1 + nparam);
 117         is_identity = isl_mat_is_scaled_identity(sub);
 118         isl_mat_free(sub);
 119
 120         return is_identity;
 121 }
 122
 123 /* Return an affine expression of the variables of the range of "morph"
 124  * in terms of the parameters and the variables of the domain on "morph".
 125  *
 126  * In order for the space manipulations to make sense, we require
 127  * that the parameters are not modified by "morph".
 128  */
 129 __isl_give isl_multi_aff *isl_morph_get_var_multi_aff(
 130         __isl_keep isl_morph *morph)
 131 {
 132         isl_space *dom, *ran, *space;
 133         isl_local_space *ls;
 134         isl_multi_aff *ma;
 135         unsigned nparam, nvar;
 136         int i;
 137         int is_identity;
 138
 139         if (!morph)
 140                 return NULL;
 141
 142         is_identity = identity_on_parameters(morph);
 143         if (is_identity < 0)
 144                 return NULL;
 145         if (!is_identity)
 146                 isl_die(isl_morph_get_ctx(morph), isl_error_invalid,
 147                         "cannot handle parameter compression", return NULL);
 148
 149         dom = isl_morph_get_dom_space(morph);
 150         ls = isl_local_space_from_space(isl_space_copy(dom));
 151         ran = isl_morph_get_ran_space(morph);
 152         space = isl_space_map_from_domain_and_range(dom, ran);
 153         ma = isl_multi_aff_zero(space);
 154
 155         nparam = isl_multi_aff_dim(ma, isl_dim_param);
 156         nvar = isl_multi_aff_dim(ma, isl_dim_out);
 157         for (i = 0; i < nvar; ++i) {
 158                 isl_val *val;
 159                 isl_vec *v;
 160                 isl_aff *aff;
 161
 162                 v = isl_mat_get_row(morph->map, 1 + nparam + i);
 163                 v = isl_vec_insert_els(v, 0, 1);
 164                 val = isl_mat_get_element_val(morph->map, 0, 0);
 165                 v = isl_vec_set_element_val(v, 0, val);
 166                 aff = isl_aff_alloc_vec(isl_local_space_copy(ls), v);
 167                 ma = isl_multi_aff_set_aff(ma, i, aff);
 168         }
 169
 170         isl_local_space_free(ls);
 171         return ma;
 172 }
 173
 174 /* Return the domain space of "morph".
 175  */
 176 __isl_give isl_space *isl_morph_get_dom_space(__isl_keep isl_morph *morph)
 177 {
 178         if (!morph)
 179                 return NULL;
 180
 181         return isl_basic_set_get_space(morph->dom);
 182 }
 183
 184 __isl_give isl_space *isl_morph_get_ran_space(__isl_keep isl_morph *morph)
 185 {
 186         if (!morph)
 187                 return NULL;
 188
 189         return isl_space_copy(morph->ran->dim);
 190 }
 191
 192 unsigned isl_morph_dom_dim(__isl_keep isl_morph *morph, enum isl_dim_type type)
 193 {
 194         if (!morph)
 195                 return 0;
 196
 197         return isl_basic_set_dim(morph->dom, type);
 198 }
 199
 200 unsigned isl_morph_ran_dim(__isl_keep isl_morph *morph, enum isl_dim_type type)
 201 {
 202         if (!morph)
 203                 return 0;
 204
 205         return isl_basic_set_dim(morph->ran, type);
 206 }
 207
 208 __isl_give isl_morph *isl_morph_remove_dom_dims(__isl_take isl_morph *morph,
 209         enum isl_dim_type type, unsigned first, unsigned n)
 210 {
 211         unsigned dom_offset;
 212
 213         if (n == 0)
 214                 return morph;
 215
 216         morph = isl_morph_cow(morph);
 217         if (!morph)
 218                 return NULL;
 219
 220         dom_offset = 1 + isl_space_offset(morph->dom->dim, type);
 221
 222         morph->dom = isl_basic_set_remove_dims(morph->dom, type, first, n);
 223
 224         morph->map = isl_mat_drop_cols(morph->map, dom_offset + first, n);
 225
 226         morph->inv = isl_mat_drop_rows(morph->inv, dom_offset + first, n);
 227
 228         if (morph->dom && morph->ran && morph->map && morph->inv)
 229                 return morph;
 230
 231         isl_morph_free(morph);
 232         return NULL;
 233 }
 234
 235 __isl_give isl_morph *isl_morph_remove_ran_dims(__isl_take isl_morph *morph,
 236         enum isl_dim_type type, unsigned first, unsigned n)
 237 {
 238         unsigned ran_offset;
 239
 240         if (n == 0)
 241                 return morph;
 242
 243         morph = isl_morph_cow(morph);
 244         if (!morph)
 245                 return NULL;
 246
 247         ran_offset = 1 + isl_space_offset(morph->ran->dim, type);
 248
 249         morph->ran = isl_basic_set_remove_dims(morph->ran, type, first, n);
 250
 251         morph->map = isl_mat_drop_rows(morph->map, ran_offset + first, n);
 252
 253         morph->inv = isl_mat_drop_cols(morph->inv, ran_offset + first, n);
 254
 255         if (morph->dom && morph->ran && morph->map && morph->inv)
 256                 return morph;
 257
 258         isl_morph_free(morph);
 259         return NULL;
 260 }
 261
 262 /* Project domain of morph onto its parameter domain.
 263  */
 264 __isl_give isl_morph *isl_morph_dom_params(__isl_take isl_morph *morph)
 265 {
 266         unsigned n;
 267
 268         morph = isl_morph_cow(morph);
 269         if (!morph)
 270                 return NULL;
 271         n = isl_basic_set_dim(morph->dom, isl_dim_set);
 272         morph = isl_morph_remove_dom_dims(morph, isl_dim_set, 0, n);
 273         if (!morph)
 274                 return NULL;
 275         morph->dom = isl_basic_set_params(morph->dom);
 276         if (morph->dom)
 277                 return morph;
 278
 279         isl_morph_free(morph);
 280         return NULL;
 281 }
 282
 283 /* Project range of morph onto its parameter domain.
 284  */
 285 __isl_give isl_morph *isl_morph_ran_params(__isl_take isl_morph *morph)
 286 {
 287         unsigned n;
 288
 289         morph = isl_morph_cow(morph);
 290         if (!morph)
 291                 return NULL;
 292         n = isl_basic_set_dim(morph->ran, isl_dim_set);
 293         morph = isl_morph_remove_ran_dims(morph, isl_dim_set, 0, n);
 294         if (!morph)
 295                 return NULL;
 296         morph->ran = isl_basic_set_params(morph->ran);
 297         if (morph->ran)
 298                 return morph;
 299
 300         isl_morph_free(morph);
 301         return NULL;
 302 }
 303
 304 void isl_morph_print_internal(__isl_take isl_morph *morph, FILE *out)
 305 {
 306         if (!morph)
 307                 return;
 308
 309         isl_basic_set_dump(morph->dom);
 310         isl_basic_set_dump(morph->ran);
 311         isl_mat_print_internal(morph->map, out, 4);
 312         isl_mat_print_internal(morph->inv, out, 4);
 313 }
 314
 315 void isl_morph_dump(__isl_take isl_morph *morph)
 316 {
 317         isl_morph_print_internal(morph, stderr);
 318 }
 319
 320 __isl_give isl_morph *isl_morph_identity(__isl_keep isl_basic_set *bset)
 321 {
 322         isl_mat *id;
 323         isl_basic_set *universe;
 324         unsigned total;
 325
 326         if (!bset)
 327                 return NULL;
 328
 329         total = isl_basic_set_total_dim(bset);
 330         id = isl_mat_identity(bset->ctx, 1 + total);
 331         universe = isl_basic_set_universe(isl_space_copy(bset->dim));
 332
 333         return isl_morph_alloc(universe, isl_basic_set_copy(universe),
 334                 id, isl_mat_copy(id));
 335 }
 336
 337 /* Create a(n identity) morphism between empty sets of the same dimension
 338  * a "bset".
 339  */
 340 __isl_give isl_morph *isl_morph_empty(__isl_keep isl_basic_set *bset)
 341 {
 342         isl_mat *id;
 343         isl_basic_set *empty;
 344         unsigned total;
 345
 346         if (!bset)
 347                 return NULL;
 348
 349         total = isl_basic_set_total_dim(bset);
 350         id = isl_mat_identity(bset->ctx, 1 + total);
 351         empty = isl_basic_set_empty(isl_space_copy(bset->dim));
 352
 353         return isl_morph_alloc(empty, isl_basic_set_copy(empty),
 354                 id, isl_mat_copy(id));
 355 }
 356
 357 /* Construct a basic set described by the "n" equalities of "bset" starting
 358  * at "first".
 359  */
 360 static __isl_give isl_basic_set *copy_equalities(__isl_keep isl_basic_set *bset,
 361         unsigned first, unsigned n)
 362 {
 363         int i, k;
 364         isl_basic_set *eq;
 365         unsigned total;
 366
 367         isl_assert(bset->ctx, bset->n_div == 0, return NULL);
 368
 369         total = isl_basic_set_total_dim(bset);
 370         eq = isl_basic_set_alloc_space(isl_space_copy(bset->dim), 0, n, 0);
 371         if (!eq)
 372                 return NULL;
 373         for (i = 0; i < n; ++i) {
 374                 k = isl_basic_set_alloc_equality(eq);
 375                 if (k < 0)
 376                         goto error;
 377                 isl_seq_cpy(eq->eq[k], bset->eq[first + i], 1 + total);
 378         }
 379
 380         return eq;
 381 error:
 382         isl_basic_set_free(eq);
 383         return NULL;
 384 }
 385
 386 /* Given a basic set, exploit the equalities in the basic set to construct
 387  * a morphism that maps the basic set to a lower-dimensional space.
 388  * Specifically, the morphism reduces the number of dimensions of type "type".
 389  *
 390  * We first select the equalities of interest, that is those that involve
 391  * variables of type "type" and no later variables.
 392  * Denote those equalities as
 393  *
 394  *              -C(p) + M x = 0
 395  *
 396  * where C(p) depends on the parameters if type == isl_dim_set and
 397  * is a constant if type == isl_dim_param.
 398  *
 399  * Use isl_mat_final_variable_compression to construct a compression
 400  *
 401  *      x = T x'
 402  *
 403  *      x' = Q x
 404  *
 405  * If T is a zero-column matrix, then the set of equality constraints
 406  * do not admit a solution.  In this case, an empty morphism is returned.
 407  *
 408  * Both matrices are extended to map the full original space to the full
 409  * compressed space.
 410  */
 411 __isl_give isl_morph *isl_basic_set_variable_compression(
 412         __isl_keep isl_basic_set *bset, enum isl_dim_type type)
 413 {
 414         unsigned otype;
 415         unsigned ntype;
 416         unsigned orest;
 417         unsigned nrest;
 418         int f_eq, n_eq;
 419         isl_space *dim;
 420         isl_mat *E, *Q, *C;
 421         isl_basic_set *dom, *ran;
 422
 423         if (!bset)
 424                 return NULL;
 425
 426         if (isl_basic_set_plain_is_empty(bset))
 427                 return isl_morph_empty(bset);
 428
 429         isl_assert(bset->ctx, bset->n_div == 0, return NULL);
 430
 431         otype = 1 + isl_space_offset(bset->dim, type);
 432         ntype = isl_basic_set_dim(bset, type);
 433         orest = otype + ntype;
 434         nrest = isl_basic_set_total_dim(bset) - (orest - 1);
 435
 436         for (f_eq = 0; f_eq < bset->n_eq; ++f_eq)
 437                 if (isl_seq_first_non_zero(bset->eq[f_eq] + orest, nrest) == -1)
 438                         break;
 439         for (n_eq = 0; f_eq + n_eq < bset->n_eq; ++n_eq)
 440                 if (isl_seq_first_non_zero(bset->eq[f_eq + n_eq] + otype, ntype) == -1)
 441                         break;
 442         if (n_eq == 0)
 443                 return isl_morph_identity(bset);
 444
 445         E = isl_mat_sub_alloc6(bset->ctx, bset->eq, f_eq, n_eq, 0, orest);
 446         C = isl_mat_final_variable_compression(E, otype - 1, &Q);
 447         if (!Q)
 448                 C = isl_mat_free(C);
 449         if (C && C->n_col == 0) {
 450                 isl_mat_free(C);
 451                 isl_mat_free(Q);
 452                 return isl_morph_empty(bset);
 453         }
 454
 455         Q = isl_mat_diagonal(Q, isl_mat_identity(bset->ctx, nrest));
 456         C = isl_mat_diagonal(C, isl_mat_identity(bset->ctx, nrest));
 457
 458         dim = isl_space_copy(bset->dim);
 459         dim = isl_space_drop_dims(dim, type, 0, ntype);
 460         dim = isl_space_add_dims(dim, type, ntype - n_eq);
 461         ran = isl_basic_set_universe(dim);
 462         dom = copy_equalities(bset, f_eq, n_eq);
 463
 464         return isl_morph_alloc(dom, ran, Q, C);
 465 }
 466
 467 /* Construct a parameter compression for "bset".
 468  * We basically just call isl_mat_parameter_compression with the right input
 469  * and then extend the resulting matrix to include the variables.
 470  *
 471  * The implementation assumes that "bset" does not have any equalities
 472  * that only involve the parameters and that isl_basic_set_gauss has
 473  * been applied to "bset".
 474  *
 475  * Let the equalities be given as
 476  *
 477  *      B(p) + A x = 0.
 478  *
 479  * We use isl_mat_parameter_compression_ext to compute the compression
 480  *
 481  *      p = T p'.
 482  */
 483 __isl_give isl_morph *isl_basic_set_parameter_compression(
 484         __isl_keep isl_basic_set *bset)
 485 {
 486         unsigned nparam;
 487         unsigned nvar;
 488         unsigned n_div;
 489         int n_eq;
 490         isl_mat *H, *B;
 491         isl_mat *map, *inv;
 492         isl_basic_set *dom, *ran;
 493
 494         if (!bset)
 495                 return NULL;
 496
 497         if (isl_basic_set_plain_is_empty(bset))
 498                 return isl_morph_empty(bset);
 499         if (bset->n_eq == 0)
 500                 return isl_morph_identity(bset);
 501
 502         n_eq = bset->n_eq;
 503         nparam = isl_basic_set_dim(bset, isl_dim_param);
 504         nvar = isl_basic_set_dim(bset, isl_dim_set);
 505         n_div = isl_basic_set_dim(bset, isl_dim_div);
 506
 507         if (isl_seq_first_non_zero(bset->eq[bset->n_eq - 1] + 1 + nparam,
 508                                     nvar + n_div) == -1)
 509                 isl_die(isl_basic_set_get_ctx(bset), isl_error_invalid,
 510                         "input not allowed to have parameter equalities",
 511                         return NULL);
 512         if (n_eq > nvar + n_div)
 513                 isl_die(isl_basic_set_get_ctx(bset), isl_error_invalid,
 514                         "input not gaussed", return NULL);
 515
 516         B = isl_mat_sub_alloc6(bset->ctx, bset->eq, 0, n_eq, 0, 1 + nparam);
 517         H = isl_mat_sub_alloc6(bset->ctx, bset->eq,
 518                                 0, n_eq, 1 + nparam, nvar + n_div);
 519         inv = isl_mat_parameter_compression_ext(B, H);
 520         inv = isl_mat_diagonal(inv, isl_mat_identity(bset->ctx, nvar));
 521         map = isl_mat_right_inverse(isl_mat_copy(inv));
 522
 523         dom = isl_basic_set_universe(isl_space_copy(bset->dim));
 524         ran = isl_basic_set_universe(isl_space_copy(bset->dim));
 525
 526         return isl_morph_alloc(dom, ran, map, inv);
 527 }
 528
 529 /* Add stride constraints to "bset" based on the inverse mapping
 530  * that was plugged in.  In particular, if morph maps x' to x,
 531  * the the constraints of the original input
 532  *
 533  *      A x' + b >= 0
 534  *
 535  * have been rewritten to
 536  *
 537  *      A inv x + b >= 0
 538  *
 539  * However, this substitution may loose information on the integrality of x',
 540  * so we need to impose that
 541  *
 542  *      inv x
 543  *
 544  * is integral.  If inv = B/d, this means that we need to impose that
 545  *
 546  *      B x = 0         mod d
 547  *
 548  * or
 549  *
 550  *      exists alpha in Z^m: B x = d alpha
 551  *
 552  * This function is similar to add_strides in isl_affine_hull.c
 553  */
 554 static __isl_give isl_basic_set *add_strides(__isl_take isl_basic_set *bset,
 555         __isl_keep isl_morph *morph)
 556 {
 557         int i, div, k;
 558         isl_int gcd;
 559
 560         if (isl_int_is_one(morph->inv->row[0][0]))
 561                 return bset;
 562
 563         isl_int_init(gcd);
 564
 565         for (i = 0; 1 + i < morph->inv->n_row; ++i) {
 566                 isl_seq_gcd(morph->inv->row[1 + i], morph->inv->n_col, &gcd);
 567                 if (isl_int_is_divisible_by(gcd, morph->inv->row[0][0]))
 568                         continue;
 569                 div = isl_basic_set_alloc_div(bset);
 570                 if (div < 0)
 571                         goto error;
 572                 isl_int_set_si(bset->div[div][0], 0);
 573                 k = isl_basic_set_alloc_equality(bset);
 574                 if (k < 0)
 575                         goto error;
 576                 isl_seq_cpy(bset->eq[k], morph->inv->row[1 + i],
 577                             morph->inv->n_col);
 578                 isl_seq_clr(bset->eq[k] + morph->inv->n_col, bset->n_div);
 579                 isl_int_set(bset->eq[k][morph->inv->n_col + div],
 580                             morph->inv->row[0][0]);
 581         }
 582
 583         isl_int_clear(gcd);
 584
 585         return bset;
 586 error:
 587         isl_int_clear(gcd);
 588         isl_basic_set_free(bset);
 589         return NULL;
 590 }
 591
 592 /* Apply the morphism to the basic set.
 593  * We basically just compute the preimage of "bset" under the inverse mapping
 594  * in morph, add in stride constraints and intersect with the range
 595  * of the morphism.
 596  */
 597 __isl_give isl_basic_set *isl_morph_basic_set(__isl_take isl_morph *morph,
 598         __isl_take isl_basic_set *bset)
 599 {
 600         isl_basic_set *res = NULL;
 601         isl_mat *mat = NULL;
 602         int i, k;
 603         int max_stride;
 604
 605         if (!morph || !bset)
 606                 goto error;
 607
 608         isl_assert(bset->ctx, isl_space_is_equal(bset->dim, morph->dom->dim),
 609                     goto error);
 610
 611         max_stride = morph->inv->n_row - 1;
 612         if (isl_int_is_one(morph->inv->row[0][0]))
 613                 max_stride = 0;
 614         res = isl_basic_set_alloc_space(isl_space_copy(morph->ran->dim),
 615                 bset->n_div + max_stride, bset->n_eq + max_stride, bset->n_ineq);
 616
 617         for (i = 0; i < bset->n_div; ++i)
 618                 if (isl_basic_set_alloc_div(res) < 0)
 619                         goto error;
 620
 621         mat = isl_mat_sub_alloc6(bset->ctx, bset->eq, 0, bset->n_eq,
 622                                         0, morph->inv->n_row);
 623         mat = isl_mat_product(mat, isl_mat_copy(morph->inv));
 624         if (!mat)
 625                 goto error;
 626         for (i = 0; i < bset->n_eq; ++i) {
 627                 k = isl_basic_set_alloc_equality(res);
 628                 if (k < 0)
 629                         goto error;
 630                 isl_seq_cpy(res->eq[k], mat->row[i], mat->n_col);
 631                 isl_seq_scale(res->eq[k] + mat->n_col, bset->eq[i] + mat->n_col,
 632                                 morph->inv->row[0][0], bset->n_div);
 633         }
 634         isl_mat_free(mat);
 635
 636         mat = isl_mat_sub_alloc6(bset->ctx, bset->ineq, 0, bset->n_ineq,
 637                                         0, morph->inv->n_row);
 638         mat = isl_mat_product(mat, isl_mat_copy(morph->inv));
 639         if (!mat)
 640                 goto error;
 641         for (i = 0; i < bset->n_ineq; ++i) {
 642                 k = isl_basic_set_alloc_inequality(res);
 643                 if (k < 0)
 644                         goto error;
 645                 isl_seq_cpy(res->ineq[k], mat->row[i], mat->n_col);
 646                 isl_seq_scale(res->ineq[k] + mat->n_col,
 647                                 bset->ineq[i] + mat->n_col,
 648                                 morph->inv->row[0][0], bset->n_div);
 649         }
 650         isl_mat_free(mat);
 651
 652         mat = isl_mat_sub_alloc6(bset->ctx, bset->div, 0, bset->n_div,
 653                                         1, morph->inv->n_row);
 654         mat = isl_mat_product(mat, isl_mat_copy(morph->inv));
 655         if (!mat)
 656                 goto error;
 657         for (i = 0; i < bset->n_div; ++i) {
 658                 isl_int_mul(res->div[i][0],
 659                                 morph->inv->row[0][0], bset->div[i][0]);
 660                 isl_seq_cpy(res->div[i] + 1, mat->row[i], mat->n_col);
 661                 isl_seq_scale(res->div[i] + 1 + mat->n_col,
 662                                 bset->div[i] + 1 + mat->n_col,
 663                                 morph->inv->row[0][0], bset->n_div);
 664         }
 665         isl_mat_free(mat);
 666
 667         res = add_strides(res, morph);
 668
 669         if (isl_basic_set_is_rational(bset))
 670                 res = isl_basic_set_set_rational(res);
 671
 672         res = isl_basic_set_simplify(res);
 673         res = isl_basic_set_finalize(res);
 674
 675         res = isl_basic_set_intersect(res, isl_basic_set_copy(morph->ran));
 676
 677         isl_morph_free(morph);
 678         isl_basic_set_free(bset);
 679         return res;
 680 error:
 681         isl_mat_free(mat);
 682         isl_morph_free(morph);
 683         isl_basic_set_free(bset);
 684         isl_basic_set_free(res);
 685         return NULL;
 686 }
 687
 688 /* Apply the morphism to the set.
 689  */
 690 __isl_give isl_set *isl_morph_set(__isl_take isl_morph *morph,
 691         __isl_take isl_set *set)
 692 {
 693         int i;
 694
 695         if (!morph || !set)
 696                 goto error;
 697
 698         isl_assert(set->ctx, isl_space_is_equal(set->dim, morph->dom->dim), goto error);
 699
 700         set = isl_set_cow(set);
 701         if (!set)
 702                 goto error;
 703
 704         isl_space_free(set->dim);
 705         set->dim = isl_space_copy(morph->ran->dim);
 706         if (!set->dim)
 707                 goto error;
 708
 709         for (i = 0; i < set->n; ++i) {
 710                 set->p[i] = isl_morph_basic_set(isl_morph_copy(morph), set->p[i]);
 711                 if (!set->p[i])
 712                         goto error;
 713         }
 714
 715         isl_morph_free(morph);
 716
 717         ISL_F_CLR(set, ISL_SET_NORMALIZED);
 718
 719         return set;
 720 error:
 721         isl_set_free(set);
 722         isl_morph_free(morph);
 723         return NULL;
 724 }
 725
 726 /* Construct a morphism that first does morph2 and then morph1.
 727  */
 728 __isl_give isl_morph *isl_morph_compose(__isl_take isl_morph *morph1,
 729         __isl_take isl_morph *morph2)
 730 {
 731         isl_mat *map, *inv;
 732         isl_basic_set *dom, *ran;
 733
 734         if (!morph1 || !morph2)
 735                 goto error;
 736
 737         map = isl_mat_product(isl_mat_copy(morph1->map), isl_mat_copy(morph2->map));
 738         inv = isl_mat_product(isl_mat_copy(morph2->inv), isl_mat_copy(morph1->inv));
 739         dom = isl_morph_basic_set(isl_morph_inverse(isl_morph_copy(morph2)),
 740                                   isl_basic_set_copy(morph1->dom));
 741         dom = isl_basic_set_intersect(dom, isl_basic_set_copy(morph2->dom));
 742         ran = isl_morph_basic_set(isl_morph_copy(morph1),
 743                                   isl_basic_set_copy(morph2->ran));
 744         ran = isl_basic_set_intersect(ran, isl_basic_set_copy(morph1->ran));
 745
 746         isl_morph_free(morph1);
 747         isl_morph_free(morph2);
 748
 749         return isl_morph_alloc(dom, ran, map, inv);
 750 error:
 751         isl_morph_free(morph1);
 752         isl_morph_free(morph2);
 753         return NULL;
 754 }
 755
 756 __isl_give isl_morph *isl_morph_inverse(__isl_take isl_morph *morph)
 757 {
 758         isl_basic_set *bset;
 759         isl_mat *mat;
 760
 761         morph = isl_morph_cow(morph);
 762         if (!morph)
 763                 return NULL;
 764
 765         bset = morph->dom;
 766         morph->dom = morph->ran;
 767         morph->ran = bset;
 768
 769         mat = morph->map;
 770         morph->map = morph->inv;
 771         morph->inv = mat;
 772
 773         return morph;
 774 }
 775
 776 /* We detect all the equalities first to avoid implicit equalities
 777  * being discovered during the computations.  In particular,
 778  * the compression on the variables could expose additional stride
 779  * constraints on the parameters.  This would result in existentially
 780  * quantified variables after applying the resulting morph, which
 781  * in turn could break invariants of the calling functions.
 782  */
 783 __isl_give isl_morph *isl_basic_set_full_compression(
 784         __isl_keep isl_basic_set *bset)
 785 {
 786         isl_morph *morph, *morph2;
 787
 788         bset = isl_basic_set_copy(bset);
 789         bset = isl_basic_set_detect_equalities(bset);
 790
 791         morph = isl_basic_set_variable_compression(bset, isl_dim_param);
 792         bset = isl_morph_basic_set(isl_morph_copy(morph), bset);
 793
 794         morph2 = isl_basic_set_parameter_compression(bset);
 795         bset = isl_morph_basic_set(isl_morph_copy(morph2), bset);
 796
 797         morph = isl_morph_compose(morph2, morph);
 798
 799         morph2 = isl_basic_set_variable_compression(bset, isl_dim_set);
 800         isl_basic_set_free(bset);
 801
 802         morph = isl_morph_compose(morph2, morph);
 803
 804         return morph;
 805 }
 806
 807 __isl_give isl_vec *isl_morph_vec(__isl_take isl_morph *morph,
 808         __isl_take isl_vec *vec)
 809 {
 810         if (!morph)
 811                 goto error;
 812
 813         vec = isl_mat_vec_product(isl_mat_copy(morph->map), vec);
 814
 815         isl_morph_free(morph);
 816         return vec;
 817 error:
 818         isl_morph_free(morph);
 819         isl_vec_free(vec);
 820         return NULL;
 821 }