1 2 /**-------------------------------------------------------------------** 3 ** CLooG ** 4 **-------------------------------------------------------------------** 5 ** constraintset.c ** 6 **-------------------------------------------------------------------** 7 ** First version: april 17th 2005 ** 8 **-------------------------------------------------------------------**/ 9 10 11/****************************************************************************** 12 * CLooG : the Chunky Loop Generator (experimental) * 13 ****************************************************************************** 14 * * 15 * Copyright (C) 2005 Cedric Bastoul * 16 * * 17 * This library is free software; you can redistribute it and/or * 18 * modify it under the terms of the GNU Lesser General Public * 19 * License as published by the Free Software Foundation; either * 20 * version 2.1 of the License, or (at your option) any later version. * 21 * * 22 * This library is distributed in the hope that it will be useful, * 23 * but WITHOUT ANY WARRANTY; without even the implied warranty of * 24 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU * 25 * Lesser General Public License for more details. * 26 * * 27 * You should have received a copy of the GNU Lesser General Public * 28 * License along with this library; if not, write to the Free Software * 29 * Foundation, Inc., 51 Franklin Street, Fifth Floor, * 30 * Boston, MA 02110-1301 USA * 31 * * 32 * CLooG, the Chunky Loop Generator * 33 * Written by Cedric Bastoul, Cedric.Bastoul@inria.fr * 34 * * 35 ******************************************************************************/ 36/* CAUTION: the english used for comments is probably the worst you ever read, 37 * please feel free to correct and improve it ! 38 */ 39 40 41# include <stdlib.h> 42# include <stdio.h> 43# include <ctype.h> 44#include <cloog/cloog.h> 45#include <cloog/matrix/constraintset.h> 46 47 48#define ALLOC(type) (type*)malloc(sizeof(type)) 49#define ALLOCN(type,n) (type*)malloc((n)*sizeof(type)) 50 51 52CloogConstraint *cloog_constraint_first(CloogConstraintSet *constraints); 53CloogConstraint *cloog_constraint_next(CloogConstraint *constraint); 54 55 56CloogConstraintSet *cloog_constraint_set_from_cloog_matrix(CloogMatrix *M) 57{ 58 return (CloogConstraintSet *)M; 59} 60 61 62void cloog_constraint_set_free(CloogConstraintSet *constraints) 63{ 64 cloog_matrix_free(&constraints->M); 65} 66 67int cloog_constraint_set_contains_level(CloogConstraintSet *constraints, 68 int level, int nb_parameters) 69{ 70 return constraints->M.NbColumns - 2 - nb_parameters >= level; 71} 72 73/* Check if the variable at position level is defined by an 74 * equality. If so, return the row number. Otherwise, return -1. 75 * 76 * If there is an equality, we can print it directly -no ambiguity-. 77 * PolyLib can give more than one equality, we use just the first one 78 * (this is a PolyLib problem, but all equalities are equivalent). 79 */ 80CloogConstraint *cloog_constraint_set_defining_equality(CloogConstraintSet *constraints, int level) 81{ 82 CloogConstraint *constraint = ALLOC(CloogConstraint); 83 int i; 84 85 constraint->set = constraints; 86 for (i = 0; i < constraints->M.NbRows; i++) 87 if (cloog_int_is_zero(constraints->M.p[i][0]) && 88 !cloog_int_is_zero(constraints->M.p[i][level])) { 89 constraint->line = &constraints->M.p[i]; 90 return constraint; 91 } 92 free(constraint); 93 return cloog_constraint_invalid(); 94} 95 96/* Check if the variable (e) at position level is defined by a 97 * pair of inequalities 98 * <a, i> + -m e + <b, p> + k1 >= 0 99 * <-a, i> + m e + <-b, p> + k2 >= 0 100 * with 0 <= k1 + k2 < m 101 * If so return the row number of the upper bound and set *lower 102 * to the row number of the lower bound. If not, return -1. 103 * 104 * If the variable at position level occurs in any other constraint, 105 * then we currently return -1. The modulo guard that we would generate 106 * would still be correct, but we would also need to generate 107 * guards corresponding to the other constraints, and this has not 108 * been implemented yet. 109 */ 110CloogConstraint *cloog_constraint_set_defining_inequalities(CloogConstraintSet *constraints, 111 int level, CloogConstraint **lower, int nb_par) 112{ 113 int i, j, k; 114 cloog_int_t m; 115 CloogMatrix *matrix = &constraints->M; 116 unsigned len = matrix->NbColumns - 2; 117 unsigned nb_iter = len - nb_par; 118 CloogConstraint *constraint; 119 120 for (i = 0; i < matrix->NbRows; i++) { 121 if (cloog_int_is_zero(matrix->p[i][level])) 122 continue; 123 if (cloog_int_is_zero(matrix->p[i][0])) 124 return cloog_constraint_invalid(); 125 if (cloog_int_is_one(matrix->p[i][level])) 126 return cloog_constraint_invalid(); 127 if (cloog_int_is_neg_one(matrix->p[i][level])) 128 return cloog_constraint_invalid(); 129 if (cloog_seq_first_non_zero(matrix->p[i]+level+1, 130 (1+nb_iter)-(level+1)) != -1) 131 return cloog_constraint_invalid(); 132 for (j = i+1; j < matrix->NbRows; ++j) { 133 if (cloog_int_is_zero(matrix->p[j][level])) 134 continue; 135 if (cloog_int_is_zero(matrix->p[j][0])) 136 return cloog_constraint_invalid(); 137 if (cloog_int_is_one(matrix->p[j][level])) 138 return cloog_constraint_invalid(); 139 if (cloog_int_is_neg_one(matrix->p[j][level])) 140 return cloog_constraint_invalid(); 141 if (cloog_seq_first_non_zero(matrix->p[j]+level+1, 142 (1+nb_iter)-(level+1)) != -1) 143 return cloog_constraint_invalid(); 144 145 cloog_int_init(m); 146 cloog_int_add(m, matrix->p[i][1+len], matrix->p[j][1+len]); 147 if (cloog_int_is_neg(m) || 148 cloog_int_abs_ge(m, matrix->p[i][level])) { 149 cloog_int_clear(m); 150 return cloog_constraint_invalid(); 151 } 152 cloog_int_clear(m); 153 154 if (!cloog_seq_is_neg(matrix->p[i]+1, matrix->p[j]+1, 155 len)) 156 return cloog_constraint_invalid(); 157 for (k = j+1; k < matrix->NbRows; ++k) 158 if (!cloog_int_is_zero(matrix->p[k][level])) 159 return cloog_constraint_invalid(); 160 *lower = ALLOC(CloogConstraint); 161 constraint = ALLOC(CloogConstraint); 162 (*lower)->set = constraints; 163 constraint->set = constraints; 164 if (cloog_int_is_pos(matrix->p[i][level])) { 165 (*lower)->line = &matrix->p[i]; 166 constraint->line = &matrix->p[j]; 167 } else { 168 (*lower)->line = &matrix->p[j]; 169 constraint->line = &matrix->p[i]; 170 } 171 return constraint; 172 } 173 } 174 return cloog_constraint_invalid(); 175} 176 177int cloog_constraint_set_total_dimension(CloogConstraintSet *constraints) 178{ 179 return constraints->M.NbColumns - 2; 180} 181 182int cloog_constraint_set_n_iterators(CloogConstraintSet *constraint, int nb_par) 183{ 184 return cloog_constraint_set_total_dimension(constraint) - nb_par; 185} 186 187int cloog_equal_total_dimension(CloogEqualities *equal) 188{ 189 return cloog_constraint_set_total_dimension(equal->constraints); 190} 191 192int cloog_constraint_total_dimension(CloogConstraint *constraint) 193{ 194 return cloog_constraint_set_total_dimension(constraint->set); 195} 196 197 198 199/****************************************************************************** 200 * Equalities spreading functions * 201 ******************************************************************************/ 202 203 204/* Equalities are stored inside a CloogMatrix data structure called "equal". 205 * This matrix has (nb_scattering + nb_iterators + 1) rows (i.e. total 206 * dimensions + 1, the "+ 1" is because a statement can be included inside an 207 * external loop without iteration domain), and (nb_scattering + nb_iterators + 208 * nb_parameters + 2) columns (all unknowns plus the scalar plus the equality 209 * type). The ith row corresponds to the equality "= 0" for the ith dimension 210 * iterator. The first column gives the equality type (0: no equality, then 211 * EQTYPE_* -see pprint.h-). At each recursion of pprint, if an equality for 212 * the current level is found, the corresponding row is updated. Then the 213 * equality if it exists is used to simplify expressions (e.g. if we have 214 * "i+1" while we know that "i=2", we simplify it in "3"). At the end of 215 * the pprint call, the corresponding row is reset to zero. 216 */ 217 218CloogEqualities *cloog_equal_alloc(int n, int nb_levels, 219 int nb_parameters) 220{ 221 int i; 222 CloogEqualities *equal = ALLOC(CloogEqualities); 223 224 equal->constraints = cloog_constraint_set_from_cloog_matrix( 225 cloog_matrix_alloc(n, nb_levels + nb_parameters + 1)); 226 equal->types = ALLOCN(int, n); 227 for (i = 0; i < n; ++i) 228 equal->types[i] = EQTYPE_NONE; 229 return equal; 230} 231 232void cloog_equal_free(CloogEqualities *equal) 233{ 234 cloog_matrix_free(&equal->constraints->M); 235 free(equal->types); 236 free(equal); 237} 238 239int cloog_equal_count(CloogEqualities *equal) 240{ 241 return equal->constraints->M.NbRows; 242} 243 244CloogConstraintSet *cloog_equal_constraints(CloogEqualities *equal) 245{ 246 return equal->constraints; 247} 248 249 250/** 251 * cloog_constraint_equal_type function : 252 * This function returns the type of the equality in the constraint (line) of 253 * (constraints) for the element (level). An equality is 'constant' iff all 254 * other factors are null except the constant one. It is a 'pure item' iff 255 * it is equal or opposite to a single variable or parameter. 256 * Otherwise it is an 'affine expression'. 257 * For instance: 258 * i = -13 is constant, i = j, j = -M are pure items, 259 * j = 2*M, i = j+1, 2*j = M are affine expressions. 260 * 261 * - constraints is the matrix of constraints, 262 * - level is the column number in equal of the element which is 'equal to', 263 ** 264 * - July 3rd 2002: first version, called pprint_equal_isconstant. 265 * - July 6th 2002: adaptation for the 3 types. 266 * - June 15th 2005: (debug) expr = domain->Constraint[line] was evaluated 267 * before checking if line != ONE_TIME_LOOP. Since 268 * ONE_TIME_LOOP is -1, an invalid read was possible. 269 * - October 19th 2005: Removal of the once-time-loop specific processing. 270 */ 271static int cloog_constraint_equal_type(CloogConstraint *constraint, int level) 272{ 273 int i, one=0 ; 274 cloog_int_t *expr; 275 276 expr = *constraint->line; 277 278 if (!cloog_int_is_one(expr[level]) && !cloog_int_is_neg_one(expr[level])) 279 return EQTYPE_EXAFFINE; 280 281 /* There is only one non null factor, and it must be +1 or -1 for 282 * iterators or parameters. 283 */ 284 for (i = 1;i <= constraint->set->M.NbColumns-2; i++) 285 if (!cloog_int_is_zero(expr[i]) && (i != level)) { 286 if ((!cloog_int_is_one(expr[i]) && !cloog_int_is_neg_one(expr[i])) || (one != 0)) 287 return EQTYPE_EXAFFINE ; 288 else 289 one = 1 ; 290 } 291 /* if the constant factor is non null, it must be alone. */ 292 if (one != 0) { 293 if (!cloog_int_is_zero(expr[constraint->set->M.NbColumns-1])) 294 return EQTYPE_EXAFFINE ; 295 } 296 else 297 return EQTYPE_CONSTANT ; 298 299 return EQTYPE_PUREITEM ; 300} 301 302 303int cloog_equal_type(CloogEqualities *equal, int level) 304{ 305 return equal->types[level-1]; 306} 307 308 309/** 310 * cloog_equal_update function: 311 * this function updates a matrix of equalities where each row corresponds to 312 * the equality "=0" of an affine expression such that the entry at column 313 * "row" (="level") is not zero. This matrix is upper-triangular, except the 314 * row number "level-1" which has to be updated for the matrix to be triangular. 315 * This function achieves the processing. 316 * - equal is the matrix to be updated, 317 * - level gives the row that has to be updated (it is actually row "level-1"), 318 * - nb_par is the number of parameters of the program. 319 ** 320 * - September 20th 2005: first version. 321 */ 322static void cloog_equal_update(CloogEqualities *equal, int level, int nb_par) 323{ int i, j ; 324 cloog_int_t gcd, factor_level, factor_outer, temp_level, temp_outer; 325 326 cloog_int_init(gcd); 327 cloog_int_init(temp_level); 328 cloog_int_init(temp_outer); 329 cloog_int_init(factor_level); 330 cloog_int_init(factor_outer); 331 332 /* For each previous level, */ 333 for (i=level-2;i>=0;i--) 334 { /* if the corresponding iterator is inside the current equality and is equal 335 * to something, 336 */ 337 if (!cloog_int_is_zero(equal->constraints->M.p[level-1][i+1]) && equal->types[i]) 338 { /* Compute the Greatest Common Divisor. */ 339 cloog_int_gcd(gcd, equal->constraints->M.p[level-1][i+1], 340 equal->constraints->M.p[i][i+1]); 341 342 /* Compute the factors to apply to each row vector element. */ 343 cloog_int_divexact(factor_level, equal->constraints->M.p[i][i+1], gcd); 344 cloog_int_divexact(factor_outer, equal->constraints->M.p[level-1][i+1], gcd); 345 346 /* Now update the row 'level'. */ 347 /* - the iterators, up to level, */ 348 for (j = 1; j <= level; j++) { 349 cloog_int_mul(temp_level, factor_level, 350 equal->constraints->M.p[level-1][j]); 351 cloog_int_mul(temp_outer, factor_outer, equal->constraints->M.p[i][j]); 352 cloog_int_sub(equal->constraints->M.p[level-1][j], temp_level, temp_outer); 353 } 354 /* - between last useful iterator (level) and the first parameter, the 355 * matrix is sparse (full of zeroes), we just do nothing there. 356 * - the parameters and the scalar. 357 */ 358 for (j = 0; j < nb_par + 1; j++) { 359 cloog_int_mul(temp_level,factor_level, 360 equal->constraints->M.p[level-1] 361 [equal->constraints->M.NbColumns-j-1]); 362 cloog_int_mul(temp_outer,factor_outer, 363 equal->constraints->M.p[i][equal->constraints->M.NbColumns-j-1]); 364 cloog_int_sub(equal->constraints->M.p[level-1] 365 [equal->constraints->M.NbColumns-j-1], 366 temp_level,temp_outer) ; 367 } 368 } 369 } 370 371 /* Normalize (divide by GCD of all elements) the updated equality. */ 372 cloog_seq_normalize(&(equal->constraints->M.p[level-1][1]), 373 equal->constraints->M.NbColumns-1); 374 375 cloog_int_clear(gcd); 376 cloog_int_clear(temp_level); 377 cloog_int_clear(temp_outer); 378 cloog_int_clear(factor_level); 379 cloog_int_clear(factor_outer); 380} 381 382 383/** 384 * cloog_equal_add function: 385 * This function updates the row (level-1) of the equality matrix (equal) with 386 * the row that corresponds to the row (line) of the matrix (matrix). 387 * - equal is the matrix of equalities, 388 * - matrix is the matrix of constraints, 389 * - level is the column number in matrix of the element which is 'equal to', 390 * - line is the line number in matrix of the constraint we want to study, 391 * - the infos structure gives the user all options on code printing and more. 392 ** 393 * - July 2nd 2002: first version. 394 * - October 19th 2005: Addition of the once-time-loop specific processing. 395 */ 396void cloog_equal_add(CloogEqualities *equal, CloogConstraintSet *constraints, 397 int level, CloogConstraint *line, int nb_par) 398{ 399 int j; 400 CloogConstraint *i = cloog_constraint_invalid(); 401 CloogMatrix *matrix = &constraints->M; 402 403 /* If we are in the case of a loop running once, this means that the equality 404 * comes from an inequality. Here we find this inequality. 405 */ 406 if (!cloog_constraint_is_valid(line)) 407 { for (i = cloog_constraint_first(constraints); 408 cloog_constraint_is_valid(i); i = cloog_constraint_next(i)) 409 if ((!cloog_int_is_zero(i->line[0][0]))&& (!cloog_int_is_zero(i->line[0][level]))) 410 { line = i ; 411 412 /* Since in once-time-loops, equalities derive from inequalities, we 413 * may have to offset the values. For instance if we have 2i>=3, the 414 * equality is in fact i=2. This may happen when the level coefficient is 415 * not 1 or -1 and the scalar value is not zero. In any other case (e.g., 416 * if the inequality is an expression including outer loop counters or 417 * parameters) the once time loop would not have been detected 418 * because of floord and ceild functions. 419 */ 420 if (cloog_int_ne_si(i->line[0][level],1) && 421 cloog_int_ne_si(i->line[0][level],-1) && 422 !cloog_int_is_zero(i->line[0][matrix->NbColumns-1])) { 423 cloog_int_t denominator; 424 425 cloog_int_init(denominator); 426 cloog_int_abs(denominator, i->line[0][level]); 427 cloog_int_fdiv_q(i->line[0][matrix->NbColumns-1], 428 i->line[0][matrix->NbColumns-1], denominator); 429 cloog_int_set_si(i->line[0][level], cloog_int_sgn(i->line[0][level])); 430 cloog_int_clear(denominator); 431 } 432 433 break ; 434 } 435 } 436 assert(cloog_constraint_is_valid(line)); 437 438 /* We update the line of equal corresponding to level: 439 * - the first element gives the equality type, 440 */ 441 equal->types[level-1] = cloog_constraint_equal_type(line, level); 442 /* - the other elements corresponding to the equality itself 443 * (the iterators up to level, then the parameters and the scalar). 444 */ 445 for (j=1;j<=level;j++) 446 cloog_int_set(equal->constraints->M.p[level-1][j], line->line[0][j]); 447 for (j = 0; j < nb_par + 1; j++) 448 cloog_int_set(equal->constraints->M.p[level-1][equal->constraints->M.NbColumns-j-1], 449 line->line[0][line->set->M.NbColumns-j-1]); 450 451 if (cloog_constraint_is_valid(i)) 452 cloog_constraint_release(line); 453 cloog_equal_update(equal, level, nb_par); 454} 455 456 457/** 458 * cloog_equal_del function : 459 * This function reset the equality corresponding to the iterator (level) 460 * in the equality matrix (equal). 461 * - July 2nd 2002: first version. 462 */ 463void cloog_equal_del(CloogEqualities *equal, int level) 464{ 465 equal->types[level-1] = EQTYPE_NONE; 466} 467 468 469 470/****************************************************************************** 471 * Processing functions * 472 ******************************************************************************/ 473 474/** 475 * Function cloog_constraint_set_normalize: 476 * This function will modify the constraint system in such a way that when 477 * there is an equality depending on the element at level 'level', there are 478 * no more (in)equalities depending on this element. For instance, try 479 * test/valilache.cloog with options -f 8 -l 9, with and without the call 480 * to this function. At a given moment, for the level L we will have 481 * 32*P=L && L>=1 (P is a lower level), this constraint system cannot be 482 * translated directly into a source code. Thus, we normalize the domain to 483 * remove L from the inequalities. In our example, this leads to 484 * 32*P=L && 32*P>=1, that can be transated to the code 485 * if (P>=1) { L=32*P ; ... }. This function solves the DaeGon Kim bug. 486 * WARNING: Remember that if there is another call to Polylib after a call to 487 * this function, we have to recall this function. 488 * -June 16th 2005: first version (adaptation from URGent June-7th-2005 by 489 * N. Vasilache). 490 * - June 21rd 2005: Adaptation for GMP. 491 * - November 4th 2005: Complete rewriting, simpler and faster. It is no more an 492 * adaptation from URGent. 493 */ 494void cloog_constraint_set_normalize(CloogConstraintSet *constraints, int level) 495{ int ref, i, j ; 496 cloog_int_t factor_i, factor_ref, temp_i, temp_ref, gcd; 497 CloogMatrix *matrix = &constraints->M; 498 499 if (matrix == NULL) 500 return ; 501 502 /* Don't "normalize" the constant term. */ 503 if (level == matrix->NbColumns-1) 504 return; 505 506 /* Let us find an equality for the current level that can be propagated. */ 507 for (ref=0;ref<matrix->NbRows;ref++) 508 if (cloog_int_is_zero(matrix->p[ref][0]) && !cloog_int_is_zero(matrix->p[ref][level])) { 509 cloog_int_init(gcd); 510 cloog_int_init(temp_i); 511 cloog_int_init(temp_ref); 512 cloog_int_init(factor_i); 513 cloog_int_init(factor_ref); 514 515 /* Row "ref" is the reference equality, now let us find a row to simplify.*/ 516 for (i=ref+1;i<matrix->NbRows;i++) 517 if (!cloog_int_is_zero(matrix->p[i][level])) { 518 /* Now let us set to 0 the "level" coefficient of row "j" using "ref". 519 * First we compute the factors to apply to each row vector element. 520 */ 521 cloog_int_gcd(gcd, matrix->p[ref][level], matrix->p[i][level]); 522 cloog_int_divexact(factor_i, matrix->p[ref][level], gcd); 523 cloog_int_divexact(factor_ref, matrix->p[i][level], gcd); 524 525 /* Maybe we are simplifying an inequality: factor_i must not be <0. */ 526 if (cloog_int_is_neg(factor_i)) { 527 cloog_int_abs(factor_i, factor_i); 528 cloog_int_neg(factor_ref, factor_ref); 529 } 530 531 /* Now update the vector. */ 532 for (j=1;j<matrix->NbColumns;j++) { 533 cloog_int_mul(temp_i, factor_i, matrix->p[i][j]); 534 cloog_int_mul(temp_ref, factor_ref, matrix->p[ref][j]); 535 cloog_int_sub(matrix->p[i][j], temp_i, temp_ref); 536 } 537 538 /* Normalize (divide by GCD of all elements) the updated vector. */ 539 cloog_seq_normalize(&(matrix->p[i][1]), matrix->NbColumns-1); 540 } 541 542 cloog_int_clear(gcd); 543 cloog_int_clear(temp_i); 544 cloog_int_clear(temp_ref); 545 cloog_int_clear(factor_i); 546 cloog_int_clear(factor_ref); 547 break ; 548 } 549} 550 551 552 553/** 554 * cloog_constraint_set_copy function: 555 * this functions builds and returns a "hard copy" (not a pointer copy) of a 556 * CloogMatrix data structure. 557 * - October 26th 2005: first version. 558 */ 559CloogConstraintSet *cloog_constraint_set_copy(CloogConstraintSet *constraints) 560{ int i, j ; 561 CloogMatrix *copy; 562 CloogMatrix *matrix = &constraints->M; 563 564 copy = cloog_matrix_alloc(matrix->NbRows, matrix->NbColumns); 565 566 for (i=0;i<matrix->NbRows;i++) 567 for (j=0;j<matrix->NbColumns;j++) 568 cloog_int_set(copy->p[i][j], matrix->p[i][j]); 569 570 return cloog_constraint_set_from_cloog_matrix(copy); 571} 572 573 574/** 575 * cloog_equal_vector_simplify function: 576 * this function simplify an affine expression with its coefficients in 577 * "vector" of length "length" thanks to an equality matrix "equal" that gives 578 * for some elements of the affine expression an equality with other elements, 579 * preferably constants. For instance, if the vector contains i+j+3 and the 580 * equality matrix gives i=n and j=2, the vector is simplified to n+3 and is 581 * returned in a new vector. 582 * - vector is the array of affine expression coefficients 583 * - equal is the matrix of equalities, 584 * - length is the vector length, 585 * - level is a level we don't want to simplify (-1 if none), 586 * - nb_par is the number of parameters of the program. 587 ** 588 * - September 20th 2005: first version. 589 * - November 2nd 2005: (debug) we are simplifying inequalities, thus we are 590 * not allowed to multiply the vector by a negative 591 * constant.Problem found after a report of Michael 592 * Classen. 593 */ 594struct cloog_vec *cloog_equal_vector_simplify(CloogEqualities *equal, cloog_int_t *vector, 595 int length, int level, int nb_par) 596{ int i, j ; 597 cloog_int_t gcd, factor_vector, factor_equal, temp_vector, temp_equal; 598 struct cloog_vec *simplified; 599 600 simplified = cloog_vec_alloc(length); 601 cloog_seq_cpy(simplified->p, vector, length); 602 603 cloog_int_init(gcd); 604 cloog_int_init(temp_vector); 605 cloog_int_init(temp_equal); 606 cloog_int_init(factor_vector); 607 cloog_int_init(factor_equal); 608 609 /* For each non-null coefficient in the vector, */ 610 for (i=length-nb_par-2;i>0;i--) 611 if (i != level) 612 { /* if the coefficient in not null, and there exists a useful equality */ 613 if ((!cloog_int_is_zero(simplified->p[i])) && equal->types[i-1]) 614 { /* Compute the Greatest Common Divisor. */ 615 cloog_int_gcd(gcd, simplified->p[i], equal->constraints->M.p[i-1][i]); 616 617 /* Compute the factors to apply to each row vector element. */ 618 cloog_int_divexact(factor_vector, equal->constraints->M.p[i-1][i], gcd); 619 cloog_int_divexact(factor_equal, simplified->p[i], gcd); 620 621 /* We are simplifying an inequality: factor_vector must not be <0. */ 622 if (cloog_int_is_neg(factor_vector)) { 623 cloog_int_abs(factor_vector, factor_vector); 624 cloog_int_neg(factor_equal, factor_equal); 625 } 626 627 /* Now update the vector. */ 628 /* - the iterators, up to the current level, */ 629 for (j=1;j<=length-nb_par-2;j++) { 630 cloog_int_mul(temp_vector, factor_vector, simplified->p[j]); 631 cloog_int_mul(temp_equal, factor_equal, equal->constraints->M.p[i-1][j]); 632 cloog_int_sub(simplified->p[j], temp_vector, temp_equal); 633 } 634 /* - between last useful iterator (i) and the first parameter, the equal 635 * matrix is sparse (full of zeroes), we just do nothing there. 636 * - the parameters and the scalar. 637 */ 638 for (j = 0; j < nb_par + 1; j++) { 639 cloog_int_mul(temp_vector, factor_vector, simplified->p[length-1-j]); 640 cloog_int_mul(temp_equal,factor_equal, 641 equal->constraints->M.p[i-1][equal->constraints->M.NbColumns-j-1]); 642 cloog_int_sub(simplified->p[length-1-j],temp_vector,temp_equal) ; 643 } 644 } 645 } 646 647 /* Normalize (divide by GCD of all elements) the updated vector. */ 648 cloog_seq_normalize(&simplified->p[1], length - 1); 649 650 cloog_int_clear(gcd); 651 cloog_int_clear(temp_vector); 652 cloog_int_clear(temp_equal); 653 cloog_int_clear(factor_vector); 654 cloog_int_clear(factor_equal); 655 656 return simplified ; 657} 658 659 660/** 661 * cloog_constraint_set_simplify function: 662 * this function simplify all constraints inside the matrix "matrix" thanks to 663 * an equality matrix "equal" that gives for some elements of the affine 664 * constraint an equality with other elements, preferably constants. 665 * For instance, if a row of the matrix contains i+j+3>=0 and the equality 666 * matrix gives i=n and j=2, the constraint is simplified to n+3>=0. The 667 * simplified constraints are returned back inside a new simplified matrix. 668 * - matrix is the set of constraints to simplify, 669 * - equal is the matrix of equalities, 670 * - level is a level we don't want to simplify (-1 if none), 671 * - nb_par is the number of parameters of the program. 672 ** 673 * - November 4th 2005: first version. 674 */ 675CloogConstraintSet *cloog_constraint_set_simplify(CloogConstraintSet *constraints, 676 CloogEqualities *equal, int level, int nb_par) 677{ int i, j, k ; 678 struct cloog_vec *vector; 679 CloogMatrix *simplified; 680 CloogMatrix *matrix = &constraints->M; 681 682 if (matrix == NULL) 683 return NULL ; 684 685 /* The simplified matrix is such that each row has been simplified thanks 686 * tho the "equal" matrix. We allocate the memory for the simplified matrix, 687 * then for each row of the original matrix, we compute the simplified 688 * vector and we copy its content into the according simplified row. 689 */ 690 simplified = cloog_matrix_alloc(matrix->NbRows, matrix->NbColumns); 691 for (i=0;i<matrix->NbRows;i++) 692 { vector = cloog_equal_vector_simplify(equal, matrix->p[i], 693 matrix->NbColumns, level, nb_par); 694 for (j=0;j<matrix->NbColumns;j++) 695 cloog_int_set(simplified->p[i][j], vector->p[j]); 696 697 cloog_vec_free(vector); 698 } 699 700 /* After simplification, it may happen that few constraints are the same, 701 * we remove them here by replacing them with 0=0 constraints. 702 */ 703 for (i=0;i<simplified->NbRows;i++) 704 for (j=i+1;j<simplified->NbRows;j++) 705 { for (k=0;k<simplified->NbColumns;k++) 706 if (cloog_int_ne(simplified->p[i][k],simplified->p[j][k])) 707 break ; 708 709 if (k == matrix->NbColumns) 710 { for (k=0;k<matrix->NbColumns;k++) 711 cloog_int_set_si(simplified->p[j][k],0); 712 } 713 } 714 715 return cloog_constraint_set_from_cloog_matrix(simplified); 716} 717 718 719/** 720 * Return clast_expr corresponding to the variable "level" (1 based) in 721 * the given constraint. 722 */ 723struct clast_expr *cloog_constraint_variable_expr(CloogConstraint *constraint, 724 int level, CloogNames *names) 725{ 726 int total_dim, nb_iter; 727 const char *name; 728 729 total_dim = cloog_constraint_total_dimension(constraint); 730 nb_iter = total_dim - names->nb_parameters; 731 732 if (level <= nb_iter) 733 name = cloog_names_name_at_level(names, level); 734 else 735 name = names->parameters[level - (nb_iter+1)] ; 736 737 return &new_clast_name(name)->expr; 738} 739 740 741/** 742 * Return true if constraint c involves variable v (zero-based). 743 */ 744int cloog_constraint_involves(CloogConstraint *constraint, int v) 745{ 746 return !cloog_int_is_zero(constraint->line[0][1+v]); 747} 748 749int cloog_constraint_is_lower_bound(CloogConstraint *constraint, int v) 750{ 751 return cloog_int_is_pos(constraint->line[0][1+v]); 752} 753 754int cloog_constraint_is_upper_bound(CloogConstraint *constraint, int v) 755{ 756 return cloog_int_is_neg(constraint->line[0][1+v]); 757} 758 759int cloog_constraint_is_equality(CloogConstraint *constraint) 760{ 761 return cloog_int_is_zero(constraint->line[0][0]); 762} 763 764void cloog_constraint_clear(CloogConstraint *constraint) 765{ 766 int k; 767 768 for (k = 1; k <= constraint->set->M.NbColumns - 2; k++) 769 cloog_int_set_si(constraint->line[0][k], 0); 770} 771 772CloogConstraintSet *cloog_constraint_set_drop_constraint( 773 CloogConstraintSet *constraints, CloogConstraint *constraint) 774{ 775 cloog_constraint_clear(constraint); 776 return constraints; 777} 778 779void cloog_constraint_coefficient_get(CloogConstraint *constraint, 780 int var, cloog_int_t *val) 781{ 782 cloog_int_set(*val, constraint->line[0][1+var]); 783} 784 785void cloog_constraint_coefficient_set(CloogConstraint *constraint, 786 int var, cloog_int_t val) 787{ 788 cloog_int_set(constraint->line[0][1+var], val); 789} 790 791void cloog_constraint_constant_get(CloogConstraint *constraint, cloog_int_t *val) 792{ 793 cloog_int_set(*val, constraint->line[0][constraint->set->M.NbColumns-1]); 794} 795 796/** 797 * Copy the coefficient of constraint c into dst in PolyLib order, 798 * i.e., first the coefficients of the variables, then the coefficients 799 * of the parameters and finally the constant. 800 */ 801void cloog_constraint_copy_coefficients(CloogConstraint *constraint, 802 cloog_int_t *dst) 803{ 804 cloog_seq_cpy(dst, constraint->line[0]+1, constraint->set->M.NbColumns-1); 805} 806 807CloogConstraint *cloog_constraint_invalid(void) 808{ 809 return NULL; 810} 811 812int cloog_constraint_is_valid(CloogConstraint *constraint) 813{ 814 return constraint != NULL; 815} 816 817 818/** 819 * Check whether there is any need for the constraint "upper" on 820 * "level" to get reduced. 821 * Yes. 822 */ 823int cloog_constraint_needs_reduction(CloogConstraint *upper, int level) 824{ 825 return 1; 826} 827 828 829/** 830 * Create a CloogConstraintSet containing enough information to perform 831 * a reduction on the upper equality (in this case lower is an invalid 832 * CloogConstraint) or the pair of inequalities upper and lower 833 * from within insert_modulo_guard. 834 * In the PolyLib backend, we return a CloogConstraintSet containting only 835 * the upper bound. The reduction will not change the stride so there 836 * will be no need to recompute the bound on the modulo expression. 837 */ 838CloogConstraintSet *cloog_constraint_set_for_reduction(CloogConstraint *upper, 839 CloogConstraint *lower) 840{ 841 CloogConstraintSet *set; 842 843 set = cloog_constraint_set_from_cloog_matrix( 844 cloog_matrix_alloc(1, upper->set->M.NbColumns)); 845 cloog_seq_cpy(set->M.p[0], upper->line[0], set->M.NbColumns); 846 return set; 847} 848 849 850/* Computes x, y and g such that g = gcd(a,b) and a*x+b*y = g */ 851static void Euclid(cloog_int_t a, cloog_int_t b, 852 cloog_int_t *x, cloog_int_t *y, cloog_int_t *g) 853{ 854 cloog_int_t c, d, e, f, tmp; 855 856 cloog_int_init(c); 857 cloog_int_init(d); 858 cloog_int_init(e); 859 cloog_int_init(f); 860 cloog_int_init(tmp); 861 cloog_int_abs(c, a); 862 cloog_int_abs(d, b); 863 cloog_int_set_si(e, 1); 864 cloog_int_set_si(f, 0); 865 while (cloog_int_is_pos(d)) { 866 cloog_int_tdiv_q(tmp, c, d); 867 cloog_int_mul(tmp, tmp, f); 868 cloog_int_sub(e, e, tmp); 869 cloog_int_tdiv_q(tmp, c, d); 870 cloog_int_mul(tmp, tmp, d); 871 cloog_int_sub(c, c, tmp); 872 cloog_int_swap(c, d); 873 cloog_int_swap(e, f); 874 } 875 cloog_int_set(*g, c); 876 if (cloog_int_is_zero(a)) 877 cloog_int_set_si(*x, 0); 878 else if (cloog_int_is_pos(a)) 879 cloog_int_set(*x, e); 880 else cloog_int_neg(*x, e); 881 if (cloog_int_is_zero(b)) 882 cloog_int_set_si(*y, 0); 883 else { 884 cloog_int_mul(tmp, a, *x); 885 cloog_int_sub(tmp, c, tmp); 886 cloog_int_divexact(*y, tmp, b); 887 } 888 cloog_int_clear(c); 889 cloog_int_clear(d); 890 cloog_int_clear(e); 891 cloog_int_clear(f); 892 cloog_int_clear(tmp); 893} 894 895/** 896 * Reduce the modulo guard expressed by "contraints" using equalities 897 * found in outer nesting levels (stored in "equal"). 898 * The modulo guard may be an equality or a pair of inequalities. 899 * In case of a pair of inequalities, "constraints" only contains the 900 * upper bound and *bound contains the bound on the 901 * corresponding modulo expression. The bound is left untouched by 902 * this function. 903 */ 904CloogConstraintSet *cloog_constraint_set_reduce(CloogConstraintSet *constraints, 905 int level, CloogEqualities *equal, int nb_par, cloog_int_t *bound) 906{ 907 int i, j, k, len, len2, nb_iter; 908 struct cloog_vec *line_vector2; 909 cloog_int_t *line, *line2, val, x, y, g; 910 911 len = constraints->M.NbColumns; 912 len2 = cloog_equal_total_dimension(equal) + 2; 913 nb_iter = len - 2 - nb_par; 914 915 cloog_int_init(val); 916 cloog_int_init(x); 917 cloog_int_init(y); 918 cloog_int_init(g); 919 920 line_vector2 = cloog_vec_alloc(len2); 921 line2 = line_vector2->p; 922 923 line = constraints->M.p[0]; 924 if (cloog_int_is_pos(line[level])) 925 cloog_seq_neg(line+1, line+1, len-1); 926 cloog_int_neg(line[level], line[level]); 927 assert(cloog_int_is_pos(line[level])); 928 929 for (i = nb_iter; i >= 1; --i) { 930 if (i == level) 931 continue; 932 cloog_int_fdiv_r(line[i], line[i], line[level]); 933 if (cloog_int_is_zero(line[i])) 934 continue; 935 936 /* Look for an earlier variable that is also a multiple of line[level] 937 * and check whether we can use the corresponding affine expression 938 * to "reduce" the modulo guard, where reduction means that we eliminate 939 * a variable, possibly at the expense of introducing other variables 940 * with smaller index. 941 */ 942 for (j = level-1; j >= 0; --j) { 943 CloogConstraint *equal_constraint; 944 if (cloog_equal_type(equal, j+1) != EQTYPE_EXAFFINE) 945 continue; 946 equal_constraint = cloog_equal_constraint(equal, j); 947 cloog_constraint_coefficient_get(equal_constraint, j, &val); 948 if (!cloog_int_is_divisible_by(val, line[level])) { 949 cloog_constraint_release(equal_constraint); 950 continue; 951 } 952 cloog_constraint_coefficient_get(equal_constraint, i-1, &val); 953 if (cloog_int_is_divisible_by(val, line[level])) { 954 cloog_constraint_release(equal_constraint); 955 continue; 956 } 957 for (k = j; k > i; --k) { 958 cloog_constraint_coefficient_get(equal_constraint, k-1, &val); 959 if (cloog_int_is_zero(val)) 960 continue; 961 if (!cloog_int_is_divisible_by(val, line[level])) 962 break; 963 } 964 if (k > i) { 965 cloog_constraint_release(equal_constraint); 966 continue; 967 } 968 cloog_constraint_coefficient_get(equal_constraint, i-1, &val); 969 Euclid(val, line[level], &x, &y, &g); 970 if (!cloog_int_is_divisible_by(val, line[i])) { 971 cloog_constraint_release(equal_constraint); 972 continue; 973 } 974 cloog_int_divexact(val, line[i], g); 975 cloog_int_neg(val, val); 976 cloog_int_mul(val, val, x); 977 cloog_int_set_si(y, 1); 978 /* Add (equal->p[j][i])^{-1} * line[i] times the equality */ 979 cloog_constraint_copy_coefficients(equal_constraint, line2+1); 980 cloog_seq_combine(line+1, y, line+1, val, line2+1, i); 981 cloog_seq_combine(line+len-nb_par-1, y, line+len-nb_par-1, 982 val, line2+len2-nb_par-1, nb_par+1); 983 cloog_constraint_release(equal_constraint); 984 break; 985 } 986 } 987 988 cloog_vec_free(line_vector2); 989 990 cloog_int_clear(val); 991 cloog_int_clear(x); 992 cloog_int_clear(y); 993 cloog_int_clear(g); 994 995 /* Make sure the line is not inverted again in the calling function. */ 996 cloog_int_neg(line[level], line[level]); 997 998 return constraints; 999} 1000 1001CloogConstraint *cloog_constraint_first(CloogConstraintSet *constraints) 1002{ 1003 CloogConstraint *c; 1004 if (constraints->M.NbRows == 0) 1005 return cloog_constraint_invalid(); 1006 c = ALLOC(CloogConstraint); 1007 c->set = constraints; 1008 c->line = &constraints->M.p[0]; 1009 return c; 1010} 1011 1012CloogConstraint *cloog_constraint_next(CloogConstraint *constraint) 1013{ 1014 constraint->line++; 1015 if (constraint->line == constraint->set->M.p + constraint->set->M.NbRows) { 1016 cloog_constraint_release(constraint); 1017 return NULL; 1018 } 1019 return constraint; 1020} 1021 1022CloogConstraint *cloog_constraint_copy(CloogConstraint *constraint) 1023{ 1024 CloogConstraint *c = ALLOC(CloogConstraint); 1025 c->set = constraint->set; 1026 c->line = constraint->line; 1027 return c; 1028} 1029 1030void cloog_constraint_release(CloogConstraint *constraint) 1031{ 1032 free(constraint); 1033} 1034 1035int cloog_constraint_set_foreach_constraint(CloogConstraintSet *constraints, 1036 int (*fn)(CloogConstraint *constraint, void *user), void *user) 1037{ 1038 CloogConstraint *c; 1039 1040 for (c = cloog_constraint_first(constraints); 1041 cloog_constraint_is_valid(c); c = cloog_constraint_next(c)) 1042 if (fn(c, user) < 0) { 1043 cloog_constraint_release(c); 1044 return -1; 1045 } 1046 1047 return 0; 1048} 1049 1050CloogConstraint *cloog_equal_constraint(CloogEqualities *equal, int j) 1051{ 1052 CloogConstraint *c = ALLOC(CloogConstraint); 1053 c->set = equal->constraints; 1054 c->line = &equal->constraints->M.p[j]; 1055 return c; 1056} 1057