2 /**-------------------------------------------------------------------**
4 **-------------------------------------------------------------------**
6 **-------------------------------------------------------------------**
7 ** First version: october 28th 2001 **
8 **-------------------------------------------------------------------**/
11 /******************************************************************************
12 * CLooG : the Chunky Loop Generator (experimental) *
13 ******************************************************************************
15 * Copyright (C) 2001-2005 Cedric Bastoul *
17 * This is free software; you can redistribute it and/or modify it under the *
18 * terms of the GNU General Public License as published by the Free Software *
19 * Foundation; either version 2 of the License, or (at your option) any later *
22 * This software is distributed in the hope that it will be useful, but *
23 * WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY *
24 * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License *
27 * You should have received a copy of the GNU General Public License along *
28 * with software; if not, write to the Free Software Foundation, Inc., *
29 * 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA *
31 * CLooG, the Chunky Loop Generator *
32 * Written by Cedric Bastoul, Cedric.Bastoul@inria.fr *
34 ******************************************************************************/
35 /* CAUTION: the english used for comments is probably the worst you ever read,
36 * please feel free to correct and improve it !
42 # include "cloog/cloog.h"
45 #ifndef USE_PPL_POWERSETS
46 # define USE_PPL_POWERSETS 1
49 /* Variables names for pretty printing. */
50 static char wild_name
[200][40];
52 static inline const char*
53 variable_output_function (ppl_dimension_type var
)
56 return wild_name
[var
+ 1];
62 error_handler (enum ppl_enum_error_code code
, const char* description
)
64 fprintf (stderr
, "PPL error code %d\n%s", code
, description
);
69 cloog_initialize (void)
71 sprintf (wild_name
[0], "1");
72 sprintf (wild_name
[1], "a");
73 sprintf (wild_name
[2], "b");
74 sprintf (wild_name
[3], "c");
75 sprintf (wild_name
[4], "d");
76 sprintf (wild_name
[5], "e");
77 sprintf (wild_name
[6], "f");
78 sprintf (wild_name
[7], "g");
79 sprintf (wild_name
[8], "h");
80 sprintf (wild_name
[9], "i");
81 sprintf (wild_name
[10], "j");
82 sprintf (wild_name
[11], "k");
83 sprintf (wild_name
[12], "l");
84 sprintf (wild_name
[13], "m");
85 sprintf (wild_name
[14], "n");
86 sprintf (wild_name
[15], "o");
87 sprintf (wild_name
[16], "p");
88 sprintf (wild_name
[17], "q");
89 sprintf (wild_name
[18], "r");
90 sprintf (wild_name
[19], "s");
91 sprintf (wild_name
[20], "t");
92 sprintf (wild_name
[21], "alpha");
93 sprintf (wild_name
[22], "beta");
94 sprintf (wild_name
[23], "gamma");
95 sprintf (wild_name
[24], "delta");
96 sprintf (wild_name
[25], "tau");
97 sprintf (wild_name
[26], "sigma");
98 sprintf (wild_name
[27], "chi");
99 sprintf (wild_name
[28], "omega");
100 sprintf (wild_name
[29], "pi");
101 sprintf (wild_name
[30], "ni");
102 sprintf (wild_name
[31], "Alpha");
103 sprintf (wild_name
[32], "Beta");
104 sprintf (wild_name
[33], "Gamma");
105 sprintf (wild_name
[34], "Delta");
106 sprintf (wild_name
[35], "Tau");
107 sprintf (wild_name
[36], "Sigma");
108 sprintf (wild_name
[37], "Chi");
109 sprintf (wild_name
[38], "Omega");
110 sprintf (wild_name
[39], "xxx");
112 if (ppl_initialize() < 0)
114 fprintf (stderr
, "Cannot initialize the Parma Polyhedra Library.\n");
118 if (ppl_restore_pre_PPL_rounding() < 0)
120 fprintf (stderr
, "Cannot restore the pre-PPL rounding mode.\n");
124 if (ppl_set_error_handler (error_handler
) < 0)
126 fprintf (stderr
, "Cannot install the custom error handler.\n");
130 if (ppl_io_set_variable_output_function (variable_output_function
) < 0)
132 fprintf (stderr
, "Cannot install the PPL custom variable output function. \n");
138 cloog_pol_copy (polyhedron pol
)
145 res
= cloog_new_pol (cloog_pol_dim (pol
), cloog_pol_nbc (pol
));
147 if (cloog_pol_nbc (pol
))
148 cloog_vector_copy (pol
->Constraint
[0], res
->Constraint
[0],
149 cloog_pol_nbc (pol
) * (cloog_pol_dim (pol
) + 2));
155 * The maximal number of rays allowed to be allocated by PolyLib. In fact since
156 * version 5.20, PolyLib automatically tune the number of rays by multiplying
157 * by 2 this number each time the maximum is reached. For unknown reasons
158 * PolyLib makes a segmentation fault if this number is too small. If this
159 * number is too small, performances will be reduced, if it is too high, memory
160 * will be saturated. Note that the option "-rays X" set this number to X.
164 /* Unused in this backend. */
166 int cloog_domain_allocated
= 0;
167 int cloog_domain_freed
= 0;
168 int cloog_domain_max
= 0;
170 /* The same for Value variables since in GMP mode they have to be freed. */
171 int cloog_value_allocated
= 0;
172 int cloog_value_freed
= 0;
173 int cloog_value_max
= 0;
177 cloog_domain_polyhedron_set (CloogDomain
* d
, polyhedra_union p
)
183 cloog_domain_set_references (CloogDomain
* d
, int i
)
189 cloog_new_domain (polyhedra_union p
)
191 CloogDomain
*domain
= (CloogDomain
*) malloc (sizeof (CloogDomain
));
193 cloog_domain_set_references (domain
, 1);
198 cloog_domain_alloc (polyhedron p
)
200 return cloog_new_domain (cloog_new_upol (p
));
203 static inline CloogDomain
*
204 cloog_check_domain_id (CloogDomain
*dom
)
209 static inline CloogDomain
*
210 cloog_check_domain (CloogDomain
*dom
)
219 cloog_vector_min_not_zero (Value
* p
, unsigned len
, int *index
, Value
* min
)
222 int i
= cloog_first_non_zero (p
, len
);
226 value_set_si (*min
, 1);
231 value_absolute (*min
, p
[i
]);
234 for (i
= i
+ 1; i
< (int) len
; i
++)
236 if (value_zero_p (p
[i
]))
239 value_absolute (x
, p
[i
]);
240 if (value_lt (x
, *min
))
242 value_assign (*min
, x
);
251 cloog_vector_gcd (Value
* p
, int len
, Value
* gcd
)
254 int i
, non_zero
, min_index
= 0;
256 q
= (Value
*) malloc (len
* sizeof (Value
));
258 for (i
= 0; i
< len
; i
++)
261 for (cp
= p
, cq
= q
, i
= 0; i
< len
; i
++, cq
++, cp
++)
262 value_absolute (*cq
, *cp
);
266 cloog_vector_min_not_zero (q
, len
, &min_index
, gcd
);
268 if (value_notone_p (*gcd
))
270 for (cq
= q
, non_zero
= 0, i
= 0; i
< len
; i
++, cq
++)
273 value_modulus (*cq
, *cq
, *gcd
);
274 non_zero
|= value_notzero_p (*cq
);
283 for (i
= 0; i
< len
; i
++)
290 cloog_matrix_combine (Value
* p1
, Value
* p2
, Value
* p3
, int x
, unsigned len
)
292 Value a1
, a2
, gcd
, b1
, b2
, n1
;
294 value_init (a1
), value_init (a2
), value_init (gcd
),
295 value_init (b1
), value_init (b2
), value_init (n1
);
297 value_assign (a1
, p1
[x
]);
298 value_assign (a2
, p2
[x
]);
300 value_absolute (b1
, a1
);
301 value_absolute (b2
, a2
);
305 value_division (a1
, a1
, gcd
);
306 value_division (a2
, a2
, gcd
);
307 value_oppose (n1
, a1
);
308 cloog_vector_combine (p1
+ 1, p2
+ 1, p3
+ 1, a2
, n1
, len
);
309 cloog_vector_normalize (p3
+ 1, len
);
311 value_clear (a1
), value_clear (a2
), value_clear (gcd
),
312 value_clear (b1
), value_clear (b2
), value_clear (n1
);
315 /* In the matrix representation an equality has a 0 in the first
316 column. When the value of the first column is 1, the row
317 represents an inequality. */
320 cloog_matrix_row_is_eq_p (CloogMatrix
*matrix
, int row
)
322 return value_zero_p (matrix
->p
[row
][0]);
325 static ppl_Constraint_t
326 cloog_build_ppl_cstr (ppl_Linear_Expression_t expr
, int ineq
)
328 ppl_Constraint_t cstr
;
333 ppl_new_Constraint (&cstr
, expr
, PPL_CONSTRAINT_TYPE_EQUAL
);
337 ppl_new_Constraint (&cstr
, expr
, PPL_CONSTRAINT_TYPE_GREATER_OR_EQUAL
);
341 /* Should not happen. */
348 /* Translates to PPL row I from MATRIX. CST is the constant part that
349 is added to the constraint. When INEQ is 1 the constraint is
350 translated as an inequality, when INEQ is 0 it is translated as an
351 equality, when INEQ has another value, the first column of the
352 matrix is read for determining the type of the constraint. */
354 static ppl_Constraint_t
355 cloog_translate_constraint (CloogMatrix
*matrix
, int i
, int cst
, int ineq
)
358 ppl_Constraint_t res
;
359 ppl_Coefficient_t coef
;
360 ppl_Linear_Expression_t expr
;
361 ppl_dimension_type dim
= matrix
->NbColumns
- 2;
365 ppl_new_Coefficient (&coef
);
366 ppl_new_Linear_Expression_with_dimension (&expr
, dim
);
368 for (j
= 1; j
< matrix
->NbColumns
- 1; j
++)
370 ppl_assign_Coefficient_from_mpz_t (coef
, matrix
->p
[i
][j
]);
371 ppl_Linear_Expression_add_to_coefficient (expr
, j
- 1, coef
);
374 value_set_si (val
, cst
);
375 value_addto (val
, matrix
->p
[i
][matrix
->NbColumns
- 1], val
);
376 ppl_assign_Coefficient_from_mpz_t (coef
, val
);
378 ppl_Linear_Expression_add_to_inhomogeneous (expr
, coef
);
379 ppl_delete_Coefficient (coef
);
381 if (ineq
!= 0 && ineq
!= 1)
382 ineq
= !cloog_matrix_row_is_eq_p (matrix
, i
);
384 res
= cloog_build_ppl_cstr (expr
, ineq
);
385 ppl_delete_Linear_Expression (expr
);
389 /* Translates to PPL the opposite of row I from MATRIX. When INEQ is
390 1 the constraint is translated as an inequality, when INEQ is 0 it
391 is translated as an equality, when INEQ has another value, the
392 first column of the matrix is read for determining the type of the
395 static ppl_Constraint_t
396 cloog_translate_oppose_constraint (CloogMatrix
*matrix
, int i
, int cst
, int ineq
)
399 ppl_Constraint_t res
;
400 ppl_Coefficient_t coef
;
401 ppl_Linear_Expression_t expr
;
402 ppl_dimension_type dim
= matrix
->NbColumns
- 2;
407 ppl_new_Coefficient (&coef
);
408 ppl_new_Linear_Expression_with_dimension (&expr
, dim
);
410 for (j
= 1; j
< matrix
->NbColumns
- 1; j
++)
412 value_oppose (val
, matrix
->p
[i
][j
]);
413 ppl_assign_Coefficient_from_mpz_t (coef
, val
);
414 ppl_Linear_Expression_add_to_coefficient (expr
, j
- 1, coef
);
417 value_oppose (val
, matrix
->p
[i
][matrix
->NbColumns
- 1]);
418 value_set_si (val1
, cst
);
419 value_addto (val
, val
, val1
);
420 ppl_assign_Coefficient_from_mpz_t (coef
, val
);
421 ppl_Linear_Expression_add_to_inhomogeneous (expr
, coef
);
422 ppl_delete_Coefficient (coef
);
424 if (ineq
!= 0 && ineq
!= 1)
425 ineq
= !cloog_matrix_row_is_eq_p (matrix
, i
);
427 res
= cloog_build_ppl_cstr (expr
, ineq
);
428 ppl_delete_Linear_Expression (expr
);
432 /* Adds to PPL the constraints from MATRIX. */
435 cloog_translate_constraint_matrix_1 (ppl_Polyhedron_t ppl
, CloogMatrix
*matrix
)
439 for (i
= 0; i
< matrix
->NbRows
; i
++)
441 ppl_Constraint_t c
= cloog_translate_constraint (matrix
, i
, 0, -1);
442 ppl_Polyhedron_add_constraint (ppl
, c
);
443 ppl_delete_Constraint (c
);
447 static ppl_Polyhedron_t
448 cloog_translate_constraint_matrix (CloogMatrix
*matrix
)
450 ppl_Polyhedron_t ppl
;
451 ppl_dimension_type dim
= matrix
->NbColumns
- 2;
453 ppl_new_C_Polyhedron_from_space_dimension (&ppl
, dim
, 0);
454 cloog_translate_constraint_matrix_1 (ppl
, matrix
);
458 /* Put the constraint matrix of polyhedron RES under Cloog's normal
459 form: Cloog expects to see
469 These two forms are equivalent but the expected form uses rightmost
470 indices for inequalities. */
473 cloog_pol_normal_form (polyhedron res
)
475 int dim
= cloog_pol_dim (res
);
476 int nrows
= cloog_pol_nbc (res
);
478 int neqs
= cloog_pol_nbeq (res
);
480 for (j
= 1; j
<= dim
; j
++)
484 for (rank
= 0; rank
< neqs
; rank
++)
485 if (j
- 1 == cloog_first_non_zero (res
->Constraint
[rank
] + 1, dim
))
487 for (i
= neqs
; i
< nrows
; i
++)
488 if (value_notzero_p (res
->Constraint
[i
][j
]))
489 cloog_matrix_combine (res
->Constraint
[i
],
490 res
->Constraint
[rank
],
491 res
->Constraint
[i
], j
, dim
+ 1);
499 cloog_translate_ppl_polyhedron_1 (ppl_Polyhedron_t pol
)
502 ppl_dimension_type dim
;
503 ppl_const_Constraint_System_t pcs
;
504 ppl_Constraint_System_const_iterator_t cit
, end
;
505 int eqs
, orig_ineqs
, ineqs
, row
, i
;
506 ppl_const_Constraint_t pc
;
508 ppl_Polyhedron_get_minimized_constraints (pol
, &pcs
);
509 ppl_new_Constraint_System_const_iterator (&cit
);
510 ppl_new_Constraint_System_const_iterator (&end
);
512 for (eqs
= 0, ineqs
= 0,
513 ppl_Constraint_System_begin (pcs
, cit
),
514 ppl_Constraint_System_end (pcs
, end
);
515 !ppl_Constraint_System_const_iterator_equal_test (cit
, end
);
516 ppl_Constraint_System_const_iterator_increment (cit
))
518 ppl_Constraint_System_const_iterator_dereference (cit
, &pc
);
519 (ppl_Constraint_type (pc
) == PPL_CONSTRAINT_TYPE_EQUAL
) ? eqs
++ : ineqs
++;
522 ppl_Polyhedron_space_dimension (pol
, &dim
);
525 if (1 || orig_ineqs
== 0)
526 res
= cloog_new_pol (dim
, eqs
+ ineqs
+ 1);
528 res
= cloog_new_pol (dim
, eqs
+ ineqs
);
531 /* Sort constraints: Cloog expects to see in matrices the equalities
532 followed by inequalities. */
535 for (ppl_Constraint_System_begin (pcs
, cit
), ppl_Constraint_System_end (pcs
, end
);
536 !ppl_Constraint_System_const_iterator_equal_test (cit
, end
);
537 ppl_Constraint_System_const_iterator_increment (cit
))
539 ppl_Coefficient_t coef
;
540 ppl_dimension_type col
;
545 ppl_new_Coefficient (&coef
);
546 ppl_Constraint_System_const_iterator_dereference (cit
, &pc
);
548 neg
= (ppl_Constraint_type (pc
) == PPL_CONSTRAINT_TYPE_LESS_THAN
549 || ppl_Constraint_type (pc
) == PPL_CONSTRAINT_TYPE_LESS_OR_EQUAL
) ? 1 : 0;
550 row
= (ppl_Constraint_type (pc
) == PPL_CONSTRAINT_TYPE_EQUAL
) ? eqs
++ : ineqs
++;
552 for (col
= 0; col
< dim
; col
++)
554 ppl_Constraint_coefficient (pc
, col
, coef
);
555 ppl_Coefficient_to_mpz_t (coef
, val
);
558 value_oppose (val
, val
);
560 value_assign (res
->Constraint
[row
][col
+ 1], val
);
563 ppl_Constraint_inhomogeneous_term (pc
, coef
);
564 ppl_Coefficient_to_mpz_t (coef
, val
);
565 value_assign (res
->Constraint
[row
][dim
+ 1], val
);
566 ppl_delete_Coefficient (coef
);
568 switch (ppl_Constraint_type (pc
))
570 case PPL_CONSTRAINT_TYPE_EQUAL
:
571 value_set_si (res
->Constraint
[row
][0], 0);
574 case PPL_CONSTRAINT_TYPE_LESS_THAN
:
575 case PPL_CONSTRAINT_TYPE_GREATER_THAN
:
576 value_decrement (res
->Constraint
[row
][dim
+ 1],
577 res
->Constraint
[row
][dim
+ 1]);
578 value_set_si (res
->Constraint
[row
][0], 1);
581 case PPL_CONSTRAINT_TYPE_LESS_OR_EQUAL
:
582 case PPL_CONSTRAINT_TYPE_GREATER_OR_EQUAL
:
583 value_set_si (res
->Constraint
[row
][0], 1);
587 fprintf (stderr
, "PPL_CONSTRAINT_TYPE_%d not implemented yet\n",
588 ppl_Constraint_type (pc
));
592 ppl_delete_Constraint_System_const_iterator (cit
);
593 ppl_delete_Constraint_System_const_iterator (end
);
595 if (cloog_pol_nbeq (res
) == 2 && cloog_pol_nbc (res
) == 2
596 && cloog_first_non_zero (res
->Constraint
[0], dim
+ 2) == (int) dim
+ 1)
598 cloog_pol_free (res
);
599 return cloog_empty_polyhedron (dim
);
602 /* Add the positivity constraint. */
603 if (1 || orig_ineqs
== 0)
606 value_set_si (res
->Constraint
[row
][0], 1);
607 for (i
= 0; i
< (int) dim
; i
++)
608 value_set_si (res
->Constraint
[row
][i
+ 1], 0);
609 value_set_si (res
->Constraint
[row
][dim
+ 1], 1);
612 /* Put the matrix of RES in normal form. */
613 cloog_pol_normal_form (res
);
615 /* If we do not sort the matrices, Cloog is a random loop
617 cloog_pol_sort_rows (res
);
623 cloog_pol_from_matrix (CloogMatrix
* m
)
626 int ncolumns
= cloog_matrix_ncolumns (m
);
627 int nrows
= cloog_matrix_nrows (m
);
631 return cloog_universe_polyhedron (ncolumns
- 2);
634 p
= cloog_translate_constraint_matrix (m
);
635 res
= cloog_translate_ppl_polyhedron_1 (p
);
636 ppl_delete_Polyhedron (p
);
637 if ((int) cloog_pol_nbc (res
) < cloog_matrix_nrows (m
))
640 cloog_pol_free (res
);
641 res
= cloog_new_pol (ncolumns
- 2, nrows
);
642 cloog_vector_copy (m
->p
[0], res
->Constraint
[0], m
->NbRows
* m
->NbColumns
);
648 cloog_domain_matrix2domain (CloogMatrix
* matrix
)
650 return cloog_domain_alloc (cloog_pol_from_matrix (matrix
));
654 cloog_translate_ppl_polyhedron (ppl_Polyhedron_t p
)
656 polyhedron res
= cloog_translate_ppl_polyhedron_1 (p
);
657 return cloog_domain_alloc (res
);
661 cloog_pol_print (FILE * file
, polyhedron pol
)
669 fprintf (file
, "<null polyhedron>\n");
673 dim
= cloog_pol_dim (pol
) + 2;
674 nbc
= cloog_pol_nbc (pol
);
675 fprintf (file
, "POLYHEDRON Dimension:%d\n", cloog_pol_dim (pol
));
677 " Constraints:%d Equations:%d\n",
678 cloog_pol_nbc (pol
), cloog_pol_nbeq (pol
));
679 fprintf (file
, "Constraints %d %d\n", nbc
, dim
);
681 for (i
= 0; i
< (int) nbc
; i
++)
683 p
= pol
->Constraint
[i
];
685 if (value_notzero_p (*p
))
686 fprintf (file
, "Inequality: [");
688 fprintf (file
, "Equality: [");
692 for (j
= 1; j
< (int) dim
; j
++)
693 value_print (file
, "%4s ", *p
++);
695 fprintf (file
, " ]\n");
699 void debug_poly (polyhedron p
)
701 cloog_pol_print (stderr
, p
);
705 debug_ppl_poly (ppl_Polyhedron_t p
)
707 debug_poly (cloog_domain_polyhedron (cloog_translate_ppl_polyhedron (p
)));
711 cloog_domain_references (CloogDomain
* d
)
713 return d
->_references
;
717 * cloog_domain_print function:
718 * This function prints the content of a CloogDomain structure (domain) into
719 * a file (foo, possibly stdout).
722 cloog_domain_print (FILE * foo
, CloogDomain
* domain
)
724 polyhedra_union upol
= cloog_domain_upol (domain
);
728 cloog_pol_print (foo
, cloog_upol_polyhedron (upol
));
729 upol
= cloog_upol_next (upol
);
732 fprintf (foo
, "Number of active references: %d\n",
733 cloog_domain_references (domain
));
737 * cloog_domain_free function:
738 * This function frees the allocated memory for a CloogDomain structure
739 * (domain). It decrements the number of active references to this structure,
740 * if there are no more references on the structure, it frees it (with the
741 * included list of polyhedra).
744 cloog_domain_free (CloogDomain
* domain
)
748 cloog_domain_set_references (domain
,
749 cloog_domain_references (domain
) - 1);
751 if (cloog_domain_references (domain
) == 0)
754 polyhedra_union upol
= cloog_domain_upol (domain
);
758 polyhedra_union next_upol
;
760 cloog_pol_free (cloog_upol_polyhedron (upol
));
761 next_upol
= cloog_upol_next (upol
);
762 cloog_upol_free (upol
);
773 * cloog_domain_copy function:
774 * This function returns a copy of a CloogDomain structure (domain). To save
775 * memory this is not a memory copy but we increment a counter of active
776 * references inside the structure, then return a pointer to that structure.
779 cloog_domain_copy (CloogDomain
* domain
)
781 cloog_domain_set_references (domain
, cloog_domain_references (domain
) + 1);
786 * cloog_domain_convex function:
787 * Given a polyhedral domain (polyhedron), this function concatenates the lists
788 * of rays and lines of the two (or more) polyhedra in the domain into one
789 * combined list, and find the set of constraints which tightly bound all of
790 * those objects. It returns the corresponding polyhedron.
793 cloog_domain_convex (CloogDomain
* domain
)
797 polyhedra_union upol
= cloog_domain_upol (domain
);
798 CloogMatrix
*m
= cloog_upol_domain2matrix (upol
);
799 ppl_Polyhedron_t p1
= cloog_translate_constraint_matrix (m
);
801 upol
= cloog_upol_next (upol
);
804 m
= cloog_upol_domain2matrix (upol
);
805 p2
= cloog_translate_constraint_matrix (m
);
806 ppl_Polyhedron_upper_bound_assign (p1
, p2
);
807 ppl_delete_Polyhedron (p2
);
809 upol
= cloog_upol_next (upol
);
812 res
= cloog_translate_ppl_polyhedron (p1
);
813 ppl_delete_Polyhedron (p1
);
819 cloog_domain_isconvex (CloogDomain
* domain
)
821 if (cloog_domain_polyhedron (domain
))
822 return !cloog_upol_next (cloog_domain_upol (domain
));
828 * cloog_domain_simple_convex:
829 * Given a list (union) of polyhedra, this function returns a "simple"
830 * convex hull of this union. In particular, the constraints of the
831 * the returned polyhedron consist of (parametric) lower and upper
832 * bounds on individual variables and constraints that appear in the
833 * original polyhedra.
835 * nb_par is the number of parameters of the domain.
838 cloog_domain_simple_convex (CloogDomain
* domain
, int nb_par
)
840 fprintf (stderr
, "cloog_domain_simple_convex (domain, nb_par = %d) is not implemented yet.\n", nb_par
);
841 fprintf (stderr
, "domain = \n");
842 cloog_domain_print (stderr
, domain
);
846 /* Returns non-zero when the constraint I in MATRIX is the positivity
847 constraint: "0 >= 0". */
850 cloog_positivity_constraint_p (CloogMatrix
*matrix
, int i
, int dim
)
854 for (j
= 1; j
< dim
; j
++)
855 if (value_notzero_p (matrix
->p
[i
][j
]))
861 /* Simplifies DOM1 in the context of DOM2. For example, DOM1 can
862 contain the following conditions: i >= 0, i <= 5, and DOM2 is a
863 loop around, i.e. the context of DOM1, that also contains the
864 conditions i >= 0, i <= 5. So instead of generating a loop like:
866 | for (i = 0; i < 6; i++)
868 | if (i >= 0 && i <= 5)
872 this function allows to detect that in the context of DOM2, the
873 constraints of DOM1 are redundant, and so the following code should
876 | for (i = 0; i < 6; i++)
880 #if USE_PPL_POWERSETS
883 cloog_domain_simplify (CloogDomain
* dom1
, CloogDomain
* dom2
)
890 CloogDomain
*res
= NULL
;
892 ppl_Pointset_Powerset_C_Polyhedron_t ps1
, ps2
;
893 ppl_dimension_type dim
= cloog_domain_dim(dom1
);
894 /* Translate dom1 into PPL powerset ps1. */
896 ppl_new_Pointset_Powerset_C_Polyhedron_from_space_dimension(&ps1
, dim
, 1);
898 for (u1
= cloog_domain_upol (dom1
); u1
; u1
= cloog_upol_next (u1
))
900 CloogMatrix
*m
= cloog_upol_domain2matrix (u1
);
901 ppl_const_Polyhedron_t ph
= cloog_translate_constraint_matrix (m
);
902 ppl_Pointset_Powerset_C_Polyhedron_add_disjunct(ps1
, ph
);
903 ppl_delete_Polyhedron(ph
);
904 cloog_matrix_free (m
);
908 /* Translate dom2 into PPL powerset ps2. */
910 ppl_new_Pointset_Powerset_C_Polyhedron_from_space_dimension(&ps2
, dim
, 1);
912 for (u2
= cloog_domain_upol (dom2
); u2
; u2
= cloog_upol_next (u2
))
914 CloogMatrix
*m
= cloog_upol_domain2matrix (u2
);
915 ppl_Polyhedron_t ph
= cloog_translate_constraint_matrix (m
);
916 ppl_Pointset_Powerset_C_Polyhedron_add_disjunct(ps2
, ph
);
917 ppl_delete_Polyhedron(ph
);
918 cloog_matrix_free (m
);
922 ppl_Pointset_Powerset_C_Polyhedron_simplify_using_context_assign(ps1
, ps2
);
924 /* Translate back simplified ps1 into res. */
925 ppl_Pointset_Powerset_C_Polyhedron_const_iterator_t i
;
926 ppl_new_Pointset_Powerset_C_Polyhedron_const_iterator(&i
);
927 ppl_Pointset_Powerset_C_Polyhedron_const_iterator_t end
;
928 ppl_new_Pointset_Powerset_C_Polyhedron_const_iterator(&end
);
929 for (ppl_Pointset_Powerset_C_Polyhedron_const_iterator_begin(ps1
, i
),
930 ppl_Pointset_Powerset_C_Polyhedron_const_iterator_end(ps1
, end
);
931 !ppl_Pointset_Powerset_C_Polyhedron_const_iterator_equal_test(i
, end
);
932 ppl_Pointset_Powerset_C_Polyhedron_const_iterator_increment(i
))
934 ppl_const_Polyhedron_t ph
;
935 ppl_Pointset_Powerset_C_Polyhedron_const_iterator_dereference(i
, &ph
);
936 res
= cloog_domain_union (res
, cloog_translate_ppl_polyhedron (ph
));
939 /* Final clean-up. */
940 ppl_delete_Pointset_Powerset_C_Polyhedron_const_iterator(i
);
941 ppl_delete_Pointset_Powerset_C_Polyhedron_const_iterator(end
);
942 ppl_delete_Pointset_Powerset_C_Polyhedron(ps1
);
943 ppl_delete_Pointset_Powerset_C_Polyhedron(ps2
);
947 #else /* !USE_PPL_POWERSETS */
950 cloog_domain_simplify (CloogDomain
* dom1
, CloogDomain
* dom2
)
952 CloogDomain
*res
= NULL
;
953 polyhedra_union u1
, u2
;
954 CloogDomain
*inter
= cloog_domain_intersection (dom1
, dom2
);
956 for (u1
= cloog_domain_upol (inter
); u1
; u1
= cloog_upol_next (u1
))
958 CloogMatrix
*m1
= cloog_upol_domain2matrix (u1
);
959 ppl_Polyhedron_t p1
= cloog_translate_constraint_matrix (m1
);
961 cloog_matrix_free (m1
);
962 for (u2
= cloog_domain_upol (dom2
); u2
; u2
= cloog_upol_next (u2
))
964 CloogMatrix
*m2
= cloog_upol_domain2matrix (u2
);
965 ppl_Polyhedron_t p2
= cloog_translate_constraint_matrix (m2
);
967 cloog_matrix_free (m2
);
968 ppl_Polyhedron_simplify_using_context_assign (p1
, p2
);
969 res
= cloog_domain_union (res
, cloog_translate_ppl_polyhedron (p1
));
970 ppl_delete_Polyhedron (p2
);
977 #endif /* !USE_PPL_POWERSETS */
979 static polyhedra_union
980 cloog_upol_copy (polyhedra_union p
)
982 polyhedra_union res
= cloog_new_upol (cloog_pol_copy (cloog_upol_polyhedron (p
)));
983 polyhedra_union upol
= res
;
985 while (cloog_upol_next (p
))
987 cloog_upol_set_next (upol
, cloog_new_upol (cloog_pol_copy (cloog_upol_polyhedron (p
))));
988 upol
= cloog_upol_next (upol
);
989 p
= cloog_upol_next (p
);
995 /* Adds to D1 the union of polyhedra from D2, with no checks of
996 redundancies between polyhedra. */
999 cloog_domain_add_domain (CloogDomain
*d1
, CloogDomain
*d2
)
1001 polyhedra_union upol
;
1009 upol
= cloog_domain_upol (d1
);
1010 while (cloog_upol_next (upol
))
1011 upol
= cloog_upol_next (upol
);
1013 cloog_upol_set_next (upol
, cloog_domain_upol (d2
));
1018 * cloog_domain_union function:
1019 * This function returns a new CloogDomain structure including a polyhedral
1020 * domain which is the union of two polyhedral domains (pol1) U (pol2) inside
1021 * two CloogDomain structures.
1024 cloog_domain_union (CloogDomain
* dom1
, CloogDomain
* dom2
)
1027 polyhedra_union head1
, head2
, tail1
, tail2
;
1028 polyhedra_union d1
, d2
;
1036 if (cloog_domain_dim (dom1
) != cloog_domain_dim (dom2
))
1038 fprintf (stderr
, "cloog_domain_union should not be called on domains of different dimensions.\n");
1044 for (d1
= cloog_domain_upol (dom1
); d1
; d1
= cloog_upol_next (d1
))
1047 ppl_Polyhedron_t ppl1
= cloog_translate_constraint_matrix (cloog_upol_domain2matrix (d1
));
1049 for (d2
= cloog_domain_upol (dom2
); d2
; d2
= cloog_upol_next (d2
))
1051 ppl_Polyhedron_t ppl2
= cloog_translate_constraint_matrix (cloog_upol_domain2matrix (d2
));
1053 if (ppl_Polyhedron_contains_Polyhedron (ppl2
, ppl1
))
1055 ppl_delete_Polyhedron (ppl2
);
1059 ppl_delete_Polyhedron (ppl2
);
1061 ppl_delete_Polyhedron (ppl1
);
1067 head1
= cloog_upol_copy (d1
);
1072 cloog_upol_set_next (tail1
, cloog_upol_copy (d1
));
1073 tail1
= cloog_upol_next (tail1
);
1080 for (d2
= cloog_domain_upol (dom2
); d2
; d2
= cloog_upol_next (d2
))
1083 ppl_Polyhedron_t ppl2
= cloog_translate_constraint_matrix (cloog_upol_domain2matrix (d2
));
1085 for (d1
= head1
; d1
; d1
= cloog_upol_next (d1
))
1087 ppl_Polyhedron_t ppl1
= cloog_translate_constraint_matrix (cloog_upol_domain2matrix (d1
));
1089 if (ppl_Polyhedron_contains_Polyhedron (ppl1
, ppl2
))
1091 ppl_delete_Polyhedron (ppl1
);
1095 ppl_delete_Polyhedron (ppl1
);
1097 ppl_delete_Polyhedron (ppl2
);
1103 head2
= cloog_upol_copy (d2
);
1108 cloog_upol_set_next (tail2
, cloog_upol_copy (d2
));
1109 tail2
= cloog_upol_next (tail2
);
1115 res
= cloog_new_domain (head2
);
1118 cloog_upol_set_next (tail1
, head2
);
1119 res
= cloog_new_domain (head1
);
1122 return cloog_check_domain (res
);
1126 * cloog_domain_intersection function:
1127 * This function returns a new CloogDomain structure including a polyhedral
1128 * domain which is the intersection of two polyhedral domains (pol1)inter(pol2)
1129 * inside two CloogDomain structures.
1132 cloog_domain_intersection (CloogDomain
* dom1
, CloogDomain
* dom2
)
1134 CloogDomain
*res
= NULL
;
1135 polyhedra_union p1
, p2
;
1136 ppl_Polyhedron_t ppl1
, ppl2
;
1138 for (p1
= cloog_domain_upol (dom1
); p1
; p1
= cloog_upol_next (p1
))
1140 ppl1
= cloog_translate_constraint_matrix (cloog_upol_domain2matrix (p1
));
1142 for (p2
= cloog_domain_upol (dom2
); p2
; p2
= cloog_upol_next (p2
))
1144 ppl2
= cloog_translate_constraint_matrix (cloog_upol_domain2matrix (p2
));
1145 ppl_Polyhedron_intersection_assign (ppl2
, ppl1
);
1146 res
= cloog_domain_union (res
, cloog_translate_ppl_polyhedron (ppl2
));
1147 ppl_delete_Polyhedron (ppl2
);
1149 ppl_delete_Polyhedron (ppl1
);
1155 /* Returns d1 minus d2. */
1158 cloog_domain_difference (CloogDomain
* d1
, CloogDomain
* d2
)
1160 CloogDomain
*res
= NULL
, *d
= d1
;
1161 polyhedra_union p1
, p2
;
1163 if (cloog_domain_isempty (d2
))
1164 return cloog_domain_copy (d1
);
1166 for (p2
= cloog_domain_upol (d2
); p2
; p2
= cloog_upol_next (p2
))
1168 CloogMatrix
*matrix
= cloog_upol_domain2matrix (p2
);
1170 for (p1
= cloog_domain_upol (d
); p1
; p1
= cloog_upol_next (p1
))
1173 CloogMatrix
*m1
= cloog_upol_domain2matrix (p1
);
1175 for (i
= 0; i
< matrix
->NbRows
; i
++)
1177 ppl_Polyhedron_t p3
;
1178 ppl_Constraint_t cstr
;
1180 /* Don't handle "0 >= 0". */
1181 if (cloog_positivity_constraint_p (matrix
, i
,
1182 cloog_domain_dim (d
) + 1))
1185 /* Add the constraint "-matrix[i] - 1 >= 0". */
1186 p3
= cloog_translate_constraint_matrix (m1
);
1187 cstr
= cloog_translate_oppose_constraint (matrix
, i
, -1, 1);
1188 ppl_Polyhedron_add_constraint (p3
, cstr
);
1189 ppl_delete_Constraint (cstr
);
1190 res
= cloog_domain_union (res
, cloog_translate_ppl_polyhedron (p3
));
1191 ppl_delete_Polyhedron (p3
);
1193 /* For an equality, add the constraint "matrix[i] - 1 >= 0". */
1194 if (cloog_matrix_row_is_eq_p (matrix
, i
))
1196 p3
= cloog_translate_constraint_matrix (m1
);
1197 cstr
= cloog_translate_constraint (matrix
, i
, -1, 1);
1198 ppl_Polyhedron_add_constraint (p3
, cstr
);
1199 ppl_delete_Constraint (cstr
);
1200 res
= cloog_domain_union (res
, cloog_translate_ppl_polyhedron (p3
));
1201 ppl_delete_Polyhedron (p3
);
1210 res
= cloog_domain_empty (cloog_domain_dim (d2
));
1219 * cloog_domain_addconstraints function :
1220 * This function adds source's polyhedron constraints to target polyhedron: for
1221 * each element of the polyhedron inside 'target' (i.e. element of the union
1222 * of polyhedra) it adds the constraints of the corresponding element in
1224 * - August 10th 2002: first version.
1225 * Nota bene for future : it is possible that source and target don't have the
1226 * same number of elements (try iftest2 without non-shared constraint
1227 * elimination in cloog_loop_separate !). This function is yet another part
1228 * of the DomainSimplify patching problem...
1231 cloog_domain_addconstraints (CloogDomain
*domain_source
, CloogDomain
*domain_target
)
1234 ppl_Polyhedron_t ppl
;
1235 polyhedra_union source
, target
, last
;
1236 int dim
= cloog_domain_dim (domain_target
);
1238 source
= cloog_domain_upol (domain_source
);
1239 target
= cloog_domain_upol (domain_target
);
1241 ppl_new_C_Polyhedron_from_space_dimension (&ppl
, dim
, 0);
1242 cloog_translate_constraint_matrix_1 (ppl
, cloog_upol_domain2matrix (target
));
1243 cloog_translate_constraint_matrix_1 (ppl
, cloog_upol_domain2matrix (source
));
1244 res
= cloog_translate_ppl_polyhedron (ppl
);
1245 ppl_delete_Polyhedron (ppl
);
1246 last
= cloog_domain_upol (res
);
1248 source
= cloog_upol_next (source
);
1249 target
= cloog_upol_next (target
);
1253 ppl_new_C_Polyhedron_from_space_dimension (&ppl
, dim
, 0);
1254 cloog_translate_constraint_matrix_1 (ppl
, cloog_upol_domain2matrix (target
));
1258 cloog_translate_constraint_matrix_1 (ppl
, cloog_upol_domain2matrix (source
));
1259 source
= cloog_upol_next (source
);
1263 (last
, cloog_domain_upol (cloog_translate_ppl_polyhedron (ppl
)));
1264 ppl_delete_Polyhedron (ppl
);
1266 last
= cloog_upol_next (last
);
1267 target
= cloog_upol_next (target
);
1273 /* Compares P1 to P2: returns 0 when the polyhedra don't overlap,
1274 returns 1 when p1 >= p2, and returns -1 when p1 < p2. The ">"
1275 relation is the "contains" relation. */
1278 cloog_domain_polyhedron_compare (CloogMatrix
*m1
, CloogMatrix
*m2
, int level
, int nb_par
, int dimension
)
1281 ppl_Polyhedron_t q1
, q2
, q3
, q4
, q5
, q
;
1282 ppl_Polyhedron_t p1
, p2
;
1284 p1
= cloog_translate_constraint_matrix (m1
);
1285 if (ppl_Polyhedron_is_empty (p1
))
1287 ppl_delete_Polyhedron (p1
);
1291 p2
= cloog_translate_constraint_matrix (m2
);
1292 if (ppl_Polyhedron_is_empty (p2
))
1294 ppl_delete_Polyhedron (p2
);
1298 ppl_new_C_Polyhedron_from_C_Polyhedron (&q1
, p1
);
1299 ppl_new_C_Polyhedron_from_C_Polyhedron (&q2
, p2
);
1301 for (i
= level
; i
< dimension
- nb_par
+ 1; i
++)
1304 ppl_Coefficient_t d
;
1305 ppl_Linear_Expression_t expr
;
1309 value_set_si (val
, 1);
1310 ppl_new_Coefficient_from_mpz_t (&d
, val
);
1311 ppl_new_Linear_Expression_with_dimension (&expr
, dimension
);
1312 ppl_Linear_Expression_add_to_coefficient (expr
, i
- 1, d
);
1313 ppl_new_Generator (&g
, expr
, PPL_GENERATOR_TYPE_LINE
, d
);
1314 ppl_Polyhedron_add_generator (q1
, g
);
1315 ppl_Polyhedron_add_generator (q2
, g
);
1316 ppl_delete_Generator (g
);
1317 ppl_delete_Linear_Expression (expr
);
1318 ppl_delete_Coefficient (d
);
1321 ppl_new_C_Polyhedron_from_C_Polyhedron (&q
, q1
);
1322 ppl_Polyhedron_intersection_assign (q
, q2
);
1323 ppl_delete_Polyhedron (q1
);
1324 ppl_delete_Polyhedron (q2
);
1326 if (ppl_Polyhedron_is_empty (q
))
1328 ppl_delete_Polyhedron (p1
);
1329 ppl_delete_Polyhedron (p2
);
1330 ppl_delete_Polyhedron (q
);
1334 ppl_new_C_Polyhedron_from_C_Polyhedron (&q1
, p1
);
1335 ppl_new_C_Polyhedron_from_C_Polyhedron (&q2
, p2
);
1336 ppl_delete_Polyhedron (p1
);
1337 ppl_delete_Polyhedron (p2
);
1339 ppl_Polyhedron_intersection_assign (q1
, q
);
1340 ppl_Polyhedron_intersection_assign (q2
, q
);
1342 m1
= cloog_upol_domain2matrix (cloog_domain_upol (cloog_translate_ppl_polyhedron (q1
)));
1343 m2
= cloog_upol_domain2matrix (cloog_domain_upol (cloog_translate_ppl_polyhedron (q2
)));
1345 ppl_new_C_Polyhedron_from_C_Polyhedron (&q4
, q
);
1346 for (i
= 0; i
< m1
->NbRows
; i
++)
1347 if (value_one_p (m1
->p
[i
][0])
1348 && value_pos_p (m1
->p
[i
][level
]))
1350 ppl_Constraint_t c
= cloog_translate_constraint (m1
, i
, 0, 1);
1351 ppl_Polyhedron_add_constraint (q4
, c
);
1352 ppl_delete_Constraint (c
);
1355 for (i
= 0; i
< m2
->NbRows
; i
++)
1356 if (value_one_p (m2
->p
[i
][0])
1357 && value_neg_p (m2
->p
[i
][level
]))
1359 ppl_Constraint_t c
= cloog_translate_constraint (m2
, i
, 0, 1);
1360 ppl_Polyhedron_add_constraint (q4
, c
);
1361 ppl_delete_Constraint (c
);
1364 if (ppl_Polyhedron_is_empty (q4
))
1366 ppl_delete_Polyhedron (q
);
1367 ppl_delete_Polyhedron (q1
);
1368 ppl_delete_Polyhedron (q2
);
1369 ppl_delete_Polyhedron (q4
);
1373 ppl_delete_Polyhedron (q4
);
1374 ppl_new_C_Polyhedron_from_C_Polyhedron (&q3
, q
);
1375 for (i
= 0; i
< m1
->NbRows
; i
++)
1377 if (value_zero_p (m1
->p
[i
][0]))
1379 if (value_zero_p (m1
->p
[i
][level
]))
1382 else if (value_neg_p (m1
->p
[i
][level
]))
1384 ppl_Constraint_t c
= cloog_translate_oppose_constraint (m1
, i
, 0, 1);
1385 ppl_Polyhedron_add_constraint (q3
, c
);
1386 ppl_delete_Constraint (c
);
1391 ppl_Constraint_t c
= cloog_translate_constraint (m1
, i
, 0, 1);
1392 ppl_Polyhedron_add_constraint (q3
, c
);
1393 ppl_delete_Constraint (c
);
1397 else if (value_neg_p (m1
->p
[i
][level
]))
1399 ppl_Constraint_t c
= cloog_translate_oppose_constraint (m1
, i
, 0, 1);
1400 ppl_Polyhedron_add_constraint (q3
, c
);
1401 ppl_delete_Constraint (c
);
1407 ppl_new_C_Polyhedron_from_C_Polyhedron (&q5
, q3
);
1408 for (j
= 0; j
< m2
->NbRows
; j
++)
1410 if (value_zero_p (m2
->p
[j
][0]))
1412 if (value_zero_p (m2
->p
[j
][level
]))
1415 else if (value_pos_p (m2
->p
[j
][level
]))
1417 ppl_Constraint_t c
= cloog_translate_oppose_constraint (m2
, j
, 0, 1);
1418 ppl_Polyhedron_add_constraint (q5
, c
);
1419 ppl_delete_Constraint (c
);
1424 ppl_Constraint_t c
= cloog_translate_constraint (m2
, j
, 0, 1);
1425 ppl_Polyhedron_add_constraint (q5
, c
);
1426 ppl_delete_Constraint (c
);
1430 else if (value_pos_p (m2
->p
[j
][level
]))
1432 ppl_Constraint_t c
= cloog_translate_oppose_constraint (m2
, j
, 0, 1);
1433 ppl_Polyhedron_add_constraint (q5
, c
);
1434 ppl_delete_Constraint (c
);
1440 if (!ppl_Polyhedron_is_empty (q5
))
1442 ppl_delete_Polyhedron (q
);
1443 ppl_delete_Polyhedron (q1
);
1444 ppl_delete_Polyhedron (q2
);
1445 ppl_delete_Polyhedron (q3
);
1446 ppl_delete_Polyhedron (q5
);
1450 /* Reinitialize Q5. */
1451 ppl_delete_Polyhedron (q5
);
1452 ppl_new_C_Polyhedron_from_C_Polyhedron (&q5
, q3
);
1455 /* Reinitialize Q3. */
1456 ppl_delete_Polyhedron (q3
);
1457 ppl_new_C_Polyhedron_from_C_Polyhedron (&q3
, q
);
1460 ppl_delete_Polyhedron (q
);
1461 ppl_delete_Polyhedron (q1
);
1462 ppl_delete_Polyhedron (q2
);
1463 ppl_delete_Polyhedron (q3
);
1464 ppl_delete_Polyhedron (q5
);
1469 * cloog_domain_sort function:
1470 * This function topologically sorts (nb_pols) polyhedra. Here (pols) is a an
1471 * array of pointers to polyhedra, (nb_pols) is the number of polyhedra,
1472 * (level) is the level to consider for partial ordering (nb_par) is the
1473 * parameter space dimension, (permut) if not NULL, is an array of (nb_pols)
1474 * integers that contains a permutation specification after call in order to
1475 * apply the topological sorting.
1479 cloog_domain_sort (CloogDomain
**doms
, unsigned nb_pols
, unsigned level
,
1480 unsigned nb_par
, int *permut
)
1483 int dim
= cloog_domain_dim (doms
[0]);
1485 for (i
= 0; i
< (int) nb_pols
; i
++)
1488 /* Note that here we do a comparison per tuple of polyhedra.
1489 PolyLib does not do this, but instead it does fewer comparisons
1490 and with a complex reasoning they infer that it some comparisons
1491 are not useful. The result is that PolyLib has wrong permutations.
1493 FIXME: In the PolyLib backend, Cloog should use this algorithm
1494 instead of PolyhedronTSort, and cloog_domain_polyhedron_compare
1495 should be implemented with a simple call to PolyhedronLTQ: these
1496 two functions produce identical answers. */
1497 for (i
= 0; i
< (int) nb_pols
; i
++)
1498 for (j
= i
+ 1; j
< (int) nb_pols
; j
++)
1500 CloogMatrix
*m1
= cloog_upol_domain2matrix (cloog_domain_upol (doms
[i
]));
1501 CloogMatrix
*m2
= cloog_upol_domain2matrix (cloog_domain_upol (doms
[j
]));
1503 if (cloog_domain_polyhedron_compare (m1
, m2
, level
, nb_par
, dim
) == 1)
1506 permut
[i
] = permut
[j
];
1513 * cloog_domain_empty function:
1514 * This function allocates the memory space for a CloogDomain structure and
1515 * sets its polyhedron to an empty polyhedron with 'dimension' dimensions.
1516 * Then it returns a pointer to the allocated space.
1517 * - June 10th 2005: first version.
1520 cloog_domain_empty (int dimension
)
1522 return cloog_domain_alloc (cloog_empty_polyhedron (dimension
));
1526 /******************************************************************************
1527 * Structure display function *
1528 ******************************************************************************/
1532 * cloog_domain_print_structure :
1533 * this function is a more human-friendly way to display the CloogDomain data
1534 * structure, it only shows the constraint system and includes an indentation
1535 * level (level) in order to work with others print_structure functions.
1536 * Written by Olivier Chorier, Luc Marchaud, Pierre Martin and Romain Tartiere.
1537 * - April 24th 2005: Initial version.
1538 * - May 26th 2005: Memory leak hunt.
1539 * - June 16th 2005: (Ced) Integration in domain.c.
1542 cloog_domain_print_structure (FILE * file
, CloogDomain
* domain
, int level
)
1545 CloogMatrix
*matrix
;
1547 /* Go to the right level. */
1548 for (i
= 0; i
< level
; i
++)
1549 fprintf (file
, "|\t");
1553 fprintf (file
, "+-- CloogDomain\n");
1555 /* Print the matrix. */
1556 matrix
= cloog_upol_domain2matrix (cloog_domain_upol (domain
));
1557 cloog_matrix_print_structure (file
, matrix
, level
);
1558 cloog_matrix_free (matrix
);
1561 for (i
= 0; i
< level
+ 1; i
++)
1562 fprintf (file
, "|\t");
1563 fprintf (file
, "\n");
1566 fprintf (file
, "+-- Null CloogDomain\n");
1572 * cloog_domain_list_print function:
1573 * This function prints the content of a CloogDomainList structure into a
1574 * file (foo, possibly stdout).
1575 * - November 6th 2001: first version.
1578 cloog_domain_list_print (FILE * foo
, CloogDomainList
* list
)
1580 while (list
!= NULL
)
1582 cloog_domain_print (foo
, cloog_domain (list
));
1583 list
= cloog_next_domain (list
);
1588 /******************************************************************************
1589 * Memory deallocation function *
1590 ******************************************************************************/
1594 * cloog_domain_list_free function:
1595 * This function frees the allocated memory for a CloogDomainList structure.
1596 * - November 6th 2001: first version.
1599 cloog_domain_list_free (CloogDomainList
* list
)
1601 CloogDomainList
*temp
;
1603 while (list
!= NULL
)
1605 temp
= cloog_next_domain (list
);
1606 cloog_domain_free (cloog_domain (list
));
1613 /******************************************************************************
1614 * Reading function *
1615 ******************************************************************************/
1619 * cloog_domain_read function:
1620 * Adaptation from the PolyLib. This function reads a matrix into a file (foo,
1621 * posibly stdin) and returns a pointer to a polyhedron containing the read
1623 * - October 18th 2001: first version.
1626 cloog_domain_read (FILE * foo
)
1628 CloogMatrix
*matrix
;
1629 CloogDomain
*domain
;
1631 matrix
= cloog_matrix_read (foo
);
1632 domain
= cloog_domain_matrix2domain (matrix
);
1633 cloog_matrix_free (matrix
);
1640 * cloog_domain_union_read function:
1641 * This function reads a union of polyhedra into a file (foo, posibly stdin) and
1642 * returns a pointer to a Polyhedron containing the read information.
1643 * - September 9th 2002: first version.
1644 * - October 29th 2005: (debug) removal of a leak counting "correction" that
1645 * was just false since ages.
1648 cloog_domain_union_read (FILE * foo
)
1650 int i
, nb_components
;
1652 CloogDomain
*domain
, *temp
, *old
;
1654 /* domain reading: nb_components (constraint matrices). */
1655 while (fgets (s
, MAX_STRING
, foo
) == 0);
1656 while ((*s
== '#' || *s
== '\n') || (sscanf (s
, " %d", &nb_components
) < 1))
1657 fgets (s
, MAX_STRING
, foo
);
1659 if (nb_components
> 0)
1660 { /* 1. first part of the polyhedra union, */
1661 domain
= cloog_domain_read (foo
);
1662 /* 2. and the nexts. */
1663 for (i
= 1; i
< nb_components
; i
++)
1664 { /* Leak counting is OK since next allocated domain is freed here. */
1665 temp
= cloog_domain_read (foo
);
1667 domain
= cloog_domain_union (temp
, old
);
1668 cloog_domain_free (temp
);
1669 cloog_domain_free (old
);
1671 return cloog_check_domain (domain
);
1679 * cloog_domain_list_read function:
1680 * This function reads a list of polyhedra into a file (foo, posibly stdin) and
1681 * returns a pointer to a CloogDomainList containing the read information.
1682 * - November 6th 2001: first version.
1685 cloog_domain_list_read (FILE * foo
)
1689 CloogDomainList
*list
, *now
, *next
;
1692 /* We read first the number of polyhedra in the list. */
1693 while (fgets (s
, MAX_STRING
, foo
) == 0);
1694 while ((*s
== '#' || *s
== '\n') || (sscanf (s
, " %d", &nb_pols
) < 1))
1695 fgets (s
, MAX_STRING
, foo
);
1697 /* Then we read the polyhedra. */
1701 list
= (CloogDomainList
*) malloc (sizeof (CloogDomainList
));
1702 cloog_set_domain (list
, cloog_domain_read (foo
));
1703 cloog_set_next_domain (list
, NULL
);
1705 for (i
= 1; i
< nb_pols
; i
++)
1707 next
= (CloogDomainList
*) malloc (sizeof (CloogDomainList
));
1708 cloog_set_domain (next
, cloog_domain_read (foo
));
1709 cloog_set_next_domain (next
, NULL
);
1710 cloog_set_next_domain (now
, next
);
1711 now
= cloog_next_domain (now
);
1718 /******************************************************************************
1719 * Processing functions *
1720 ******************************************************************************/
1723 * cloog_domain_isempty function:
1724 * This function returns 1 if the polyhedron given as input is empty, 0
1726 * - October 28th 2001: first version.
1729 cloog_domain_isempty (CloogDomain
* domain
)
1731 if (cloog_domain_polyhedron (domain
) == NULL
)
1734 if (cloog_upol_next (cloog_domain_upol (domain
)))
1737 return (cloog_domain_dim (domain
) < cloog_domain_nbeq (domain
)) ? 1 : 0;
1741 * cloog_domain_project function:
1742 * From Quillere's LoopGen 0.4. This function returns the projection of
1743 * (domain) on the (level) first dimensions (i.e. outer loops). It returns a
1744 * pointer to the projected Polyhedron.
1745 * - nb_par is the number of parameters.
1747 * - October 27th 2001: first version.
1748 * - June 21rd 2005: Adaptation for GMP (based on S. Verdoolaege's version of
1752 cloog_domain_project (CloogDomain
* domain
, int level
, int nb_par
)
1754 CloogDomain
*res
= NULL
;
1755 polyhedra_union upol
= cloog_domain_upol (domain
);
1756 int i
, diff
= cloog_domain_dim (domain
) - level
- nb_par
;
1758 ppl_dimension_type
*to_remove
;
1762 fprintf (stderr
, "cloog_domain_project should not be called with"
1763 "cloog_domain_dim (domain) < level + nb_par");
1768 return cloog_domain_copy (domain
);
1771 to_remove
= (ppl_dimension_type
*) malloc (n
* sizeof (ppl_dimension_type
));
1773 for (i
= 0; i
< n
; i
++)
1774 to_remove
[i
] = i
+ level
;
1778 ppl_Polyhedron_t ppl
= cloog_translate_constraint_matrix (cloog_upol_domain2matrix (upol
));
1780 ppl_Polyhedron_remove_space_dimensions (ppl
, to_remove
, n
);
1781 res
= cloog_domain_add_domain (res
, cloog_translate_ppl_polyhedron (ppl
));
1782 ppl_delete_Polyhedron (ppl
);
1783 upol
= cloog_upol_next (upol
);
1790 * cloog_domain_extend function:
1791 * From Quillere's LoopGen 0.4. This function returns the (domain) given as
1792 * input with (dim)+(nb_par) dimensions. The new dimensions are added before
1793 * the (nb_par) parameters. This function does not free (domain), and returns
1795 * - nb_par is the number of parameters.
1797 * - October 27th 2001: first version.
1798 * - June 21rd 2005: Adaptation for GMP (based on S. Verdoolaege's version of
1802 cloog_domain_extend (CloogDomain
* domain
, int dim
, int nb_par
)
1804 CloogDomain
*res
= NULL
;
1805 polyhedra_union upol
= cloog_domain_upol (domain
);
1806 int i
, k
, n
, diff
= dim
+ nb_par
- cloog_domain_dim (domain
);
1807 ppl_dimension_type
*map
;
1808 ppl_dimension_type to_add
= diff
;
1811 return cloog_domain_copy (domain
);
1814 map
= (ppl_dimension_type
*) malloc (n
* sizeof (ppl_dimension_type
));
1816 k
= cloog_domain_dim (domain
) - nb_par
;
1817 for (i
= 0; i
< k
; i
++)
1826 map
[i
] = i
- nb_par
;
1830 ppl_Polyhedron_t ppl
= cloog_translate_constraint_matrix (cloog_upol_domain2matrix (upol
));
1832 ppl_Polyhedron_add_space_dimensions_and_embed (ppl
, to_add
);
1833 ppl_Polyhedron_map_space_dimensions (ppl
, map
, n
);
1834 res
= cloog_domain_add_domain (res
, cloog_translate_ppl_polyhedron (ppl
));
1835 ppl_delete_Polyhedron (ppl
);
1836 upol
= cloog_upol_next (upol
);
1843 * cloog_domain_never_integral function:
1844 * For us, an equality like 3*i -4 = 0 is always false since 4%3 != 0. This
1845 * function returns a boolean set to 1 if there is this kind of 'never true'
1846 * constraint inside a polyhedron, 0 otherwise.
1847 * - domain is the polyhedron to check,
1849 * - November 28th 2001: first version.
1850 * - June 26th 2003: for iterators, more 'never true' constraints are found
1851 * (compare cholesky2 and vivien with a previous version),
1852 * checking for the parameters created (compare using vivien).
1853 * - June 28th 2003: Previously in loop.c and called
1854 * cloog_loop_simplify_nevertrue, now here !
1855 * - June 21rd 2005: Adaptation for GMP (based on S. Verdoolaege's version of
1857 * - October 14th 2005: Complete rewriting, not faster but code quite shorter.
1860 cloog_domain_never_integral (CloogDomain
* domain
)
1862 int i
, dimension
, nbc
;
1866 if ((domain
== NULL
) || (cloog_domain_polyhedron (domain
) == NULL
))
1870 value_init_c (modulo
);
1871 p
= cloog_domain_polyhedron (domain
);
1872 dimension
= cloog_domain_dim (domain
) + 2;
1873 nbc
= cloog_domain_nbconstraints (domain
);
1875 /* For each constraint... */
1876 for (i
= 0; i
< nbc
; i
++)
1877 { /* If we have an equality and the scalar part is not zero... */
1878 if (value_zero_p (p
->Constraint
[i
][0]) &&
1879 value_notzero_p (p
->Constraint
[i
][dimension
- 1]))
1880 { /* Then we check whether the scalar can be divided by the gcd of the
1881 * unknown vector (including iterators and parameters) or not. If not,
1882 * there is no integer point in the polyhedron and we return 1.
1884 cloog_vector_gcd (&(p
->Constraint
[i
][1]), dimension
- 2, &gcd
);
1885 value_modulus (modulo
,
1886 p
->Constraint
[i
][dimension
- 1],
1889 if (value_notzero_p (modulo
))
1891 value_clear_c (gcd
);
1892 value_clear_c (modulo
);
1898 value_clear_c (gcd
);
1899 value_clear_c (modulo
);
1904 cloog_vector_matrix_smallest_column (Value
* a
, int n
, int p
, int q
)
1906 int res
, numero
= 0, k
, allzero
;
1912 c
= a
+ q
- 1 + p
* (n
- 1);
1914 value_absolute (comp
, *c
);
1916 for (k
= n
; k
> q
; k
--, c
-= p
, value_absolute (comp
, *c
))
1918 if (value_notzero_p (comp
))
1922 value_assign (minus
, comp
);
1926 else if (value_ge (minus
, comp
))
1928 value_assign (minus
, comp
);
1938 value_absolute (comp
, *c
);
1939 if ((value_notzero_p (comp
)) && (value_ge (minus
, comp
)))
1945 value_clear (minus
);
1952 cloog_vector_matrix_swap_rows (Value
* a
, int i
, int j
, int p
)
1955 Value
*c1
= a
+ (i
- 1) * p
;
1956 Value
*c2
= a
+ (j
- 1) * p
;
1961 for (k
= 1; k
<= p
; k
++)
1963 value_assign (s
, *c1
);
1964 value_assign (*c1
, *c2
);
1965 value_assign (*c2
, s
);
1974 cloog_vector_matrix_swap_columns (Value
* a
, int i
, int j
, int n
, int p
)
1984 for (k
= 1; k
<= n
; k
++)
1986 value_assign (s
, *c1
);
1987 value_assign (*c1
, *c2
);
1988 value_assign (*c2
, s
);
1997 cloog_vector_matrix_oppose_line (Value
* a
, int i
, int p
)
2000 Value
*c
= a
+ (i
- 1) * p
;
2002 for (k
= 1; k
<= p
; k
++, c
++)
2003 value_oppose (*c
, *c
);
2007 cloog_vector_matrix_oppose_column (Value
* a
, int i
, int n
, int p
)
2010 Value
*c
= a
+ (i
- 1);
2012 for (k
= 1; k
<= n
; k
++, c
+= p
)
2013 value_oppose (*c
, *c
);
2017 cloog_vector_matrix_combine_line (Value
* a
, int i
, int j
,
2021 Value
*c1
= a
+ (i
- 1) * p
;
2022 Value
*c2
= a
+ (j
- 1) * p
;
2024 for (k
= 1; k
<= p
; k
++, c1
++, c2
++)
2025 value_addmul (*c1
, x
, *c2
);
2029 cloog_vector_matrix_combine_column (Value
* a
, int i
, int j
, Value x
, int n
, int p
)
2032 Value
*c1
= a
+ (i
- 1);
2033 Value
*c2
= a
+ (j
- 1);
2035 for (k
= 1; k
<= n
; k
++, c1
= c1
+ p
, c2
= c2
+ p
)
2036 value_addmul (*c1
, x
, *c2
);
2040 cloog_matrix_hermite_combine (Value
*a
, Value
*b
, Value c
, Value
*d
, int k
, int n
,
2041 int p
, int q
, Value pivot
, Value x
, Value x_inv
)
2046 value_division (x
, c
, pivot
);
2047 value_modulus (tmp
, c
, pivot
);
2049 if (value_neg_p (tmp
))
2050 value_decrement (x
, x
);
2054 value_oppose (x_inv
, x
);
2055 cloog_vector_matrix_combine_line (a
, k
, q
, x_inv
, p
);
2056 cloog_vector_matrix_combine_column (b
, q
, k
, x
, n
, n
);
2057 cloog_vector_matrix_combine_line (d
, k
, q
, x_inv
, n
);
2061 cloog_matrix_hermite_oppose (Value
*a
, Value
*b
, Value
*d
, int n
, int p
, int q
, Value pivot
)
2063 cloog_vector_matrix_oppose_line (a
, q
, p
);
2064 cloog_vector_matrix_oppose_line (d
, q
, n
);
2065 cloog_vector_matrix_oppose_column (b
, q
, n
, n
);
2066 value_oppose (pivot
, pivot
);
2070 cloog_matrix_hermite_1 (Value
* a
, Value
* b
, Value
* d
, int n
, int p
, int q
)
2073 Value x
, pivot
, x_inv
;
2083 for (i
= cloog_vector_matrix_smallest_column (a
, n
, p
, q
); i
!= 0;
2084 i
= cloog_vector_matrix_smallest_column (a
, n
, p
, q
))
2088 cloog_vector_matrix_swap_rows (a
, i
, q
, p
);
2089 cloog_vector_matrix_swap_columns (b
, i
, q
, n
, n
);
2090 cloog_vector_matrix_swap_rows (d
, i
, q
, n
);
2093 c
= a
+ (q
- 1) * (p
+ 1);
2094 value_assign (pivot
, *c
);
2095 if (value_neg_p (pivot
))
2096 cloog_matrix_hermite_oppose (a
, b
, d
, n
, p
, q
, pivot
);
2098 for (c
+= p
, k
= q
+ 1; k
<= n
; k
++, c
+= p
)
2099 if (value_notzero_p (*c
))
2100 cloog_matrix_hermite_combine (a
, b
, *c
, d
, k
, n
, p
, q
, pivot
, x
, x_inv
);
2104 value_assign (pivot
, *(c
+ (q
- 1) * p
));
2105 if (value_neg_p (pivot
))
2106 cloog_matrix_hermite_oppose (a
, b
, d
, n
, p
, q
, pivot
);
2108 if (value_notzero_p (pivot
))
2109 for (k
= 1; k
< q
; k
++, c
+= p
)
2110 if (value_notzero_p (*c
))
2111 cloog_matrix_hermite_combine (a
, b
, *c
, d
, k
, n
, p
, q
, pivot
, x
, x_inv
);
2113 cloog_matrix_hermite_1 (a
, b
, d
, n
, p
, q
+ 1);
2115 value_clear (pivot
);
2117 value_clear (x_inv
);
2121 static CloogMatrix
*
2122 cloog_matrix_add_null_row (CloogMatrix
* M
)
2125 CloogMatrix
*res
= cloog_matrix_alloc (M
->NbRows
+ 1, M
->NbColumns
);
2127 for (i
= 0; i
< M
->NbRows
; i
++)
2128 for (j
= 0; j
< M
->NbColumns
; j
++)
2129 value_assign (res
->p
[i
][j
], M
->p
[i
][j
]);
2131 for (j
= 0; j
< M
->NbColumns
; j
++)
2132 value_set_si (res
->p
[i
][j
], 0);
2137 static CloogMatrix
*
2138 cloog_matrix_transpose (CloogMatrix
* m
)
2141 CloogMatrix
*res
= cloog_matrix_alloc (m
->NbColumns
, m
->NbRows
);
2143 for (i
= 0; i
< (int) m
->NbRows
; i
++)
2144 for (j
= 0; j
< (int) m
->NbColumns
; j
++)
2145 value_assign (res
->p
[j
][i
], m
->p
[i
][j
]);
2151 cloog_vector_matrix_vectorify (CloogMatrix
* A
)
2154 Value
*res
= (Value
*) malloc (sizeof (Value
) * A
->NbRows
* A
->NbColumns
);
2156 for (i
= 0; i
< A
->NbRows
* A
->NbColumns
; i
++)
2157 value_init (res
[i
]);
2159 for (i
= 0; i
< A
->NbRows
; i
++)
2160 for (j
= 0; j
< A
->NbColumns
; j
++)
2161 value_assign (res
[i
* A
->NbColumns
+ j
], A
->p
[i
][j
]);
2166 static CloogMatrix
*
2167 cloog_vector_matrix_matrixify (Value
* A
, int NbRows
, int NbCols
)
2170 CloogMatrix
*res
= cloog_matrix_alloc (NbRows
, NbCols
);
2172 for (i
= 0; i
< NbRows
; i
++)
2173 for (j
= 0; j
< NbCols
; j
++)
2174 value_assign (res
->p
[i
][j
], A
[i
* NbCols
+ j
]);
2180 cloog_vector_matrix_identity (int n
)
2183 Value
*res
= (Value
*) malloc (sizeof (Value
) * n
* n
);
2186 for (i
= 0; i
< n
* n
; i
++)
2187 value_init (res
[i
]);
2189 for (i
= 1; i
<= n
; i
++)
2190 for (j
= 1; j
<= n
; j
++, ptr
++)
2192 value_set_si (*ptr
, 1);
2194 value_set_si (*ptr
, 0);
2200 cloog_matrix_hermite (CloogMatrix
* a
, CloogMatrix
** H
, CloogMatrix
** U
)
2203 CloogMatrix
*tu
, *transpose
= cloog_matrix_transpose (a
);
2204 Value
*va
= cloog_vector_matrix_vectorify (transpose
);
2205 Value
*vid
= cloog_vector_matrix_identity (a
->NbColumns
);
2206 Value
*vid_inv
= cloog_vector_matrix_identity (a
->NbColumns
);
2208 cloog_matrix_free (transpose
);
2210 cloog_matrix_hermite_1 (va
, vid
, vid_inv
,
2211 a
->NbColumns
, a
->NbRows
, 1);
2213 transpose
= cloog_vector_matrix_matrixify (va
, a
->NbColumns
, a
->NbRows
);
2214 H
[0] = cloog_matrix_transpose (transpose
);
2215 cloog_matrix_free (transpose
);
2217 tu
= cloog_vector_matrix_matrixify (vid
, a
->NbColumns
, a
->NbColumns
);
2218 U
[0] = cloog_matrix_transpose (tu
);
2219 cloog_matrix_free (tu
);
2221 for (i
= 0; i
< a
->NbRows
* a
->NbColumns
; i
++)
2222 value_clear (va
[i
]);
2224 for (i
= 0; i
< a
->NbColumns
* a
->NbColumns
; i
++)
2225 value_clear (vid
[i
]);
2227 for (i
= 0; i
< a
->NbColumns
* a
->NbColumns
; i
++)
2228 value_clear (vid_inv
[i
]);
2236 cloog_exchange_rows (CloogMatrix
* m
, int row1
, int row2
)
2240 for (i
= 0; i
< (int) m
->NbColumns
; i
++)
2241 value_swap (m
->p
[row1
][i
], m
->p
[row2
][i
]);
2244 static CloogMatrix
*
2245 cloog_dio_copy_submatrix (CloogMatrix
* matrix
)
2248 CloogMatrix
* copy
= cloog_matrix_alloc (matrix
->NbRows
- 1, matrix
->NbColumns
- 1) ;
2250 for (i
= 0; i
< copy
->NbRows
; i
++)
2251 for (j
= 0; j
< copy
->NbColumns
; j
++)
2252 value_assign(copy
->p
[i
][j
], matrix
->p
[i
][j
]) ;
2258 cloog_dio_rearrange_redundant_rows (CloogMatrix
* M
)
2260 int i
, end
, row
, rank
= 1;
2261 CloogMatrix
*A
, *L
, *H
, *U
;
2263 L
= cloog_dio_copy_submatrix (M
);
2267 A
= cloog_matrix_alloc (1, L
->NbColumns
);
2268 for (i
= 0; i
< L
->NbColumns
; i
++)
2269 value_assign (A
->p
[0][i
], L
->p
[0][i
]);
2271 end
= L
->NbRows
- 1;
2277 CloogMatrix
*temp
= cloog_matrix_add_null_row (A
);
2279 for (i
= 0; i
< A
->NbColumns
; i
++)
2280 value_assign (temp
->p
[row
][i
], L
->p
[row
][i
]);
2282 cloog_matrix_hermite (temp
, &H
, &U
);
2283 cloog_matrix_free (U
);
2286 for (i
= 0; i
< n
; i
++)
2287 if (value_zero_p (H
->p
[i
][i
]))
2290 cloog_matrix_free (H
);
2294 /* This is a redundant row: put it at the end. */
2295 cloog_exchange_rows (M
, row
, end
);
2302 cloog_matrix_free (A
);
2303 A
= cloog_matrix_copy (temp
);
2304 cloog_matrix_free (temp
);
2307 if (rank
== (end
>= L
->NbColumns
? L
->NbColumns
: end
)
2312 cloog_matrix_free (A
);
2315 cloog_matrix_free (L
);
2320 cloog_matrix_inverse_pivot (CloogMatrix
*m
, int low
, int up
, int i
, int j
,
2321 Value m1
, Value m2
, Value x
, Value piv
)
2325 for (c
= low
; c
< up
; ++c
)
2327 value_multiply (m1
, piv
, m
->p
[j
][c
]);
2328 value_multiply (m2
, x
, m
->p
[i
][c
]);
2329 value_subtract (m
->p
[j
][c
], m1
, m2
);
2334 cloog_matrix_inverse_den (CloogMatrix
*Mat
, CloogMatrix
*MatInv
, int k
, Value
*x
)
2338 Value
*res
= (Value
*) malloc (k
* sizeof (Value
));
2343 value_set_si (*x
, 1);
2345 for (j
= 0; j
< cloog_matrix_nrows (Mat
); ++j
)
2347 value_init (res
[j
]);
2348 value_assign (res
[j
], Mat
->p
[j
][j
]);
2349 cloog_vector_gcd (&MatInv
->p
[j
][0], k
, &gcd
);
2350 Gcd (gcd
, res
[j
], &gcd
);
2352 if (value_neg_p (res
[j
]))
2353 value_oppose (gcd
, gcd
);
2355 if (value_notone_p (gcd
))
2357 for (c
= 0; c
< k
; c
++)
2358 value_division (MatInv
->p
[j
][c
], MatInv
->p
[j
][c
], gcd
);
2359 value_division (res
[j
], res
[j
], gcd
);
2362 Gcd (*x
, res
[j
], &gcd
);
2363 value_division (m1
, res
[j
], gcd
);
2364 value_multiply (*x
, *x
, m1
);
2373 cloog_matrix_inverse_div (CloogMatrix
*MatInv
, Value
*den
, Value m1
, Value x
)
2377 if (value_notone_p (x
))
2378 for (j
= 0; j
< MatInv
->NbRows
; ++j
)
2380 for (c
= 0; c
< MatInv
->NbColumns
; c
++)
2382 value_division (m1
, x
, den
[j
]);
2383 value_multiply (MatInv
->p
[j
][c
], MatInv
->p
[j
][c
], m1
);
2389 cloog_matrix_inverse_triang (CloogMatrix
*Mat
, CloogMatrix
*MatInv
, Value
*m1
, Value
*m2
)
2399 for (i
= 0; i
< cloog_matrix_ncolumns (Mat
); ++i
)
2401 if (value_zero_p (Mat
->p
[i
][i
]))
2403 for (j
= i
; j
< cloog_matrix_nrows (Mat
); ++j
)
2404 if (value_notzero_p (Mat
->p
[j
][i
]))
2407 if (j
== cloog_matrix_nrows (Mat
))
2410 cloog_matrix_exchange_rows (Mat
, i
, j
);
2411 cloog_matrix_exchange_rows (MatInv
, i
, j
);
2414 for (j
= 0; j
< cloog_matrix_nrows (Mat
) ; ++j
)
2417 || !value_notzero_p (Mat
->p
[j
][i
]))
2420 value_assign (x
, Mat
->p
[j
][i
]);
2421 value_assign (piv
, Mat
->p
[i
][i
]);
2423 if (value_notone_p (gcd
))
2425 value_division (x
, x
, gcd
);
2426 value_division (piv
, piv
, gcd
);
2429 cloog_matrix_inverse_pivot (Mat
, ((j
> i
) ? i
: 0),
2430 cloog_matrix_nrows (Mat
),
2431 i
, j
, *m1
, *m2
, x
, piv
);
2432 cloog_matrix_inverse_pivot (MatInv
, 0, cloog_matrix_nrows (Mat
),
2433 i
, j
, *m1
, *m2
, x
, piv
);
2435 cloog_vector_gcd (&MatInv
->p
[j
][0], (unsigned) cloog_matrix_nrows (Mat
), m1
);
2436 cloog_vector_gcd (&Mat
->p
[j
][0], (unsigned) cloog_matrix_nrows (Mat
), m2
);
2437 Gcd (*m1
, *m2
, &gcd
);
2438 if (value_notone_p (gcd
))
2442 for (c
= 0; c
< cloog_matrix_nrows (Mat
); ++c
)
2444 value_division (Mat
->p
[j
][c
], Mat
->p
[j
][c
], gcd
);
2445 value_division (MatInv
->p
[j
][c
], MatInv
->p
[j
][c
], gcd
);
2461 cloog_matrix_inverse (CloogMatrix
* Mat
, CloogMatrix
* MatInv
)
2469 if (Mat
->NbRows
!= Mat
->NbColumns
)
2476 cloog_vector_set (MatInv
->p
[0], 0, cloog_matrix_nrows (Mat
) * cloog_matrix_nrows (Mat
));
2478 for (i
= 0; i
< cloog_matrix_nrows (Mat
); ++i
)
2479 value_set_si (MatInv
->p
[i
][i
], 1);
2481 if (!cloog_matrix_inverse_triang (Mat
, MatInv
, &m1
, &m2
))
2484 den
= cloog_matrix_inverse_den (Mat
, MatInv
, cloog_matrix_nrows (Mat
), &x
);
2485 cloog_matrix_inverse_div (MatInv
, den
, m1
, x
);
2489 for (j
= 0; j
< cloog_matrix_nrows (Mat
); ++j
)
2490 value_clear (den
[j
]);
2502 cloog_dio_copy (CloogMatrix
*m1
, CloogMatrix
*m2
)
2505 int n
= m2
->NbRows
- 1;
2506 int m
= m2
->NbColumns
- 1;
2508 for (i
= 0; i
< n
; i
++)
2509 for (j
= 0; j
< m
; j
++)
2510 value_assign (m1
->p
[i
][j
], m2
->p
[i
][j
]);
2514 cloog_dio_negate_coeffs (CloogMatrix
*m
)
2517 int n
= m
->NbRows
- 1;
2518 int k
= m
->NbColumns
- 1;
2519 Value
*res
= (Value
*) malloc (sizeof (Value
) * (n
));
2521 for (i
= 0; i
< n
; i
++)
2523 value_init (res
[i
]);
2524 value_oppose (res
[i
], m
->p
[i
][k
]);
2531 cloog_dio_get_first_diagonal_zero (CloogMatrix
*m
)
2533 int i
, n
= m
->NbRows
<= m
->NbColumns
? m
->NbRows
: m
->NbColumns
;
2536 for (i
= 0; i
< n
; i
++)
2537 if (value_notzero_p (m
->p
[i
][i
]))
2546 cloog_dio_reduce_diagonal (CloogMatrix
*m
, Value
*coeffs
, int nbc
, int rank
)
2549 Value
*res
= (Value
*) malloc (sizeof (Value
) * nbc
);
2555 for (i
= 0; i
< nbc
; i
++)
2556 value_init (res
[i
]);
2558 for (i
= 0; i
< rank
; i
++)
2560 value_set_si (sum
, 0);
2561 for (j
= 0; j
< i
; j
++)
2562 value_addmul (sum
, res
[j
], m
->p
[i
][j
]);
2564 value_subtract (tmp
, coeffs
[i
], sum
);
2565 value_modulus (tmp
, tmp
, m
->p
[i
][i
]);
2566 if (value_notzero_p (tmp
))
2568 for (i
= 0; i
< nbc
; i
++)
2569 value_clear (res
[i
]);
2578 value_subtract (tmp
, coeffs
[i
], sum
);
2579 value_division (res
[i
], tmp
, m
->p
[i
][i
]);
2582 for (i
= rank
; i
< m
->NbColumns
; i
++)
2583 value_set_si (res
[i
], 0);
2589 cloog_dio_check_coeffs (CloogMatrix
*m
, Value
*T
, Value
*coeffs
, int nbr
, int nbc
, int rank
)
2597 for (i
= rank
; i
< m
->NbRows
; i
++)
2599 value_set_si (sum
, 0);
2600 for (j
= 0; j
< m
->NbColumns
; j
++)
2601 value_addmul (sum
, T
[j
], m
->p
[i
][j
]);
2603 if (value_ne (sum
, coeffs
[i
]))
2605 for (i
= 0; i
< nbr
; i
++)
2606 value_clear (coeffs
[i
]);
2607 for (i
= 0; i
< nbc
; i
++)
2622 static CloogMatrix
*
2623 cloog_dio_init_U (CloogMatrix
*u
, int n
, int rank
)
2629 return cloog_matrix_alloc (0, 0);
2631 res
= cloog_matrix_alloc (n
, n
- rank
);
2633 for (i
= 0; i
< res
->NbRows
; i
++)
2634 for (j
= 0; j
< res
->NbColumns
; j
++)
2635 value_assign (res
->p
[i
][j
], u
->p
[i
][j
+ rank
]);
2641 cloog_dio_init_X (CloogMatrix
*u
, Value
*t
, int n
)
2644 Vector
*res
= Vector_Alloc (n
);
2648 for (i
= 0; i
< u
->NbRows
; i
++)
2650 value_set_si (sum
, 0);
2652 for (j
= 0; j
< u
->NbColumns
; j
++)
2653 value_addmul (sum
, u
->p
[i
][j
], t
[j
]);
2655 value_assign (res
->p
[i
], sum
);
2663 cloog_solve_diophantine (CloogMatrix
* m
, CloogMatrix
** u
, Vector
** x
)
2665 int i
, nbr
, nbc
, rank
;
2666 CloogMatrix
*A
, *temp
, *hermi
, *unimod
, *unimodinv
;
2670 A
= cloog_matrix_copy (m
);
2671 nbr
= A
->NbRows
- 1;
2673 cloog_dio_rearrange_redundant_rows (A
);
2675 temp
= cloog_matrix_alloc (A
->NbRows
- 1, A
->NbColumns
- 1);
2676 cloog_dio_copy (temp
, A
);
2677 coeffs
= cloog_dio_negate_coeffs (A
);
2678 cloog_matrix_free (A
);
2680 cloog_matrix_hermite (temp
, &hermi
, &unimod
);
2681 rank
= cloog_dio_get_first_diagonal_zero (hermi
);
2682 nbc
= temp
->NbColumns
;
2684 T
= cloog_dio_reduce_diagonal (hermi
, coeffs
, nbc
, rank
);
2685 unimodinv
= cloog_matrix_alloc (unimod
->NbRows
, unimod
->NbColumns
);
2688 || !cloog_dio_check_coeffs (hermi
, T
, coeffs
, nbr
, nbc
, rank
)
2689 || !cloog_matrix_inverse (unimod
, unimodinv
))
2691 *u
= cloog_matrix_alloc (0, 0);
2692 *x
= Vector_Alloc (0);
2694 for (i
= 0; i
< nbr
; i
++)
2695 value_clear (coeffs
[i
]);
2701 *u
= cloog_dio_init_U (unimodinv
, hermi
->NbColumns
, rank
);
2702 *x
= cloog_dio_init_X (unimodinv
, T
, m
->NbColumns
- 1);
2704 cloog_matrix_free (unimod
);
2705 cloog_matrix_free (unimodinv
);
2706 cloog_matrix_free (hermi
);
2707 cloog_matrix_free (temp
);
2709 for (i
= 0; i
< nbr
; i
++)
2710 value_clear (coeffs
[i
]);
2712 for (i
= 0; i
< nbc
; i
++)
2721 * cloog_domain_stride function:
2722 * This function finds the stride imposed to unknown with the column number
2723 * 'strided_level' in order to be integral. For instance, if we have a
2724 * constraint like -i - 2j + 2k = 0, and we consider k, then k can be integral
2725 * only if (i + 2j)%2 = 0. Then only if i%2 = 0. Then k imposes a stride 2 to
2726 * the unknown i. The function returns the imposed stride in a parameter field.
2727 * - domain is the set of constraint we have to consider,
2728 * - strided_level is the column number of the unknown for which a stride have
2730 * - looking_level is the column number of the unknown that impose a stride to
2731 * the first unknown.
2732 * - stride is the stride that is returned back as a function parameter.
2733 * - offset is the value of the constant c if the condition is of the shape
2734 * (i + c)%s = 0, s being the stride.
2736 * - June 28th 2003: first version.
2737 * - July 14th 2003: can now look for multiple striding constraints and returns
2738 * the GCD of the strides and the common offset.
2739 * - June 21rd 2005: Adaptation for GMP (based on S. Verdoolaege's version of
2743 cloog_domain_stride (CloogDomain
*domain
, int strided_level
, int nb_par
, Value
*stride
, Value
*offset
)
2746 int n_col
, n_row
, rank
;
2750 polyhedron p
= cloog_domain_polyhedron (domain
);
2751 int dimension
= cloog_domain_dim (domain
);
2752 int nbeq
= cloog_domain_nbeq (domain
);
2754 /* Look at all equalities involving strided_level and the inner
2755 * iterators. We can ignore the outer iterators and the parameters
2756 * here because the lower bound on strided_level is assumed to
2759 n_col
= (1 + dimension
- nb_par
) - strided_level
;
2760 for (i
= 0, n_row
= 0; i
< nbeq
; i
++)
2761 if (cloog_first_non_zero
2762 (p
->Constraint
[i
] + strided_level
, n_col
) != -1)
2765 M
= cloog_matrix_alloc (n_row
+ 1, n_col
+ 1);
2766 for (i
= 0, n_row
= 0; i
< nbeq
; i
++)
2768 if (cloog_first_non_zero
2769 (p
->Constraint
[i
] + strided_level
, n_col
) == -1)
2771 cloog_vector_copy (p
->Constraint
[i
] + strided_level
,
2772 M
->p
[n_row
], n_col
);
2773 value_assign (M
->p
[n_row
][n_col
],
2774 p
->Constraint
[i
][1 + dimension
]);
2777 value_set_si (M
->p
[n_row
][n_col
], 1);
2779 /* Then look at the general solution to the above equalities. */
2780 rank
= cloog_solve_diophantine (M
, &U
, &V
);
2781 cloog_matrix_free (M
);
2785 /* There is no solution, so the body of this loop will
2786 * never execute. We just leave the constraints alone here so
2787 * that they will ensure the body will not be executed.
2788 * We should probably propagate this information up so that
2789 * the loop can be removed entirely.
2791 value_set_si (*offset
, 0);
2792 value_set_si (*stride
, 1);
2796 /* Compute the gcd of the coefficients defining strided_level. */
2797 cloog_vector_gcd (U
->p
[0], U
->NbColumns
, stride
);
2798 value_oppose (*offset
, V
->p
[0]);
2799 value_pmodulus (*offset
, *offset
, *stride
);
2802 cloog_matrix_free (U
);
2809 * cloog_domain_integral_lowerbound function:
2810 * This function returns 1 if the lower bound of an iterator (such as its
2811 * column rank in the constraint set 'domain' is 'level') is integral,
2812 * 0 otherwise. If the lower bound is actually integral, the function fills
2813 * the 'lower' field with the lower bound value.
2814 * - June 29th 2003: first version.
2815 * - June 21rd 2005: Adaptation for GMP (based on S. Verdoolaege's version of
2819 cloog_domain_integral_lowerbound (CloogDomain
*domain
, int level
, Value
*lower
)
2821 int i
, first_lower
= 1, dimension
, lower_constraint
= -1, nbc
;
2822 Value iterator
, constant
, tmp
;
2823 polyhedron p
= cloog_domain_polyhedron (domain
);
2824 dimension
= cloog_domain_dim (domain
);
2825 nbc
= cloog_domain_nbconstraints (domain
);
2827 /* We want one and only one lower bound (e.g. no equality, no maximum
2830 for (i
= 0; i
< nbc
; i
++)
2831 if (value_zero_p (p
->Constraint
[i
][0])
2832 && value_notzero_p (p
->Constraint
[i
][level
]))
2835 for (i
= 0; i
< nbc
; i
++)
2836 if (value_pos_p (p
->Constraint
[i
][level
]))
2841 lower_constraint
= i
;
2850 /* We want an integral lower bound: no other non-zero entry except the
2851 * iterator coefficient and the constant.
2853 for (i
= 1; i
< level
; i
++)
2855 (p
->Constraint
[lower_constraint
][i
]))
2858 for (i
= level
+ 1; i
<= dimension
; i
++)
2860 (p
->Constraint
[lower_constraint
][i
]))
2865 value_init_c (iterator
);
2866 value_init_c (constant
);
2869 /* If all is passed, then find the lower bound and return 1. */
2870 value_assign (iterator
,
2871 p
->Constraint
[lower_constraint
][level
]);
2872 value_oppose (constant
,
2873 p
->Constraint
[lower_constraint
][dimension
+ 1]);
2875 value_modulus (tmp
, constant
, iterator
);
2876 value_division (*lower
, constant
, iterator
);
2878 if (!(value_zero_p (tmp
) || value_neg_p (constant
)))
2879 value_increment (*lower
, *lower
);
2881 value_clear_c (iterator
);
2882 value_clear_c (constant
);
2883 value_clear_c (tmp
);
2889 * cloog_domain_lowerbound_update function:
2890 * This function updates the integral lower bound of an iterator (such as its
2891 * column rank in the constraint set 'domain' is 'level') into 'lower'.
2892 * - Jun 29th 2003: first version.
2893 * - June 21rd 2005: Adaptation for GMP (based on S. Verdoolaege's version of
2897 cloog_domain_lowerbound_update (CloogDomain
*domain
, int level
, Value lower
)
2900 int nbc
= cloog_domain_nbconstraints (domain
);
2901 int dim
= cloog_domain_dim (domain
);
2902 polyhedron p
= cloog_domain_polyhedron (domain
);
2904 /* There is only one lower bound, the first one is the good one. */
2905 for (i
= 0; i
< nbc
; i
++)
2906 if (value_pos_p (p
->Constraint
[i
][level
]))
2908 value_set_si (p
->Constraint
[i
][level
], 1);
2909 value_oppose (p
->Constraint
[i
][dim
+ 1], lower
);
2916 * cloog_domain_lazy_equal function:
2917 * This function returns 1 if the domains given as input are the same, 0 if it
2918 * is unable to decide. This function makes an entry-to-entry comparison between
2919 * the constraint systems, if all the entries are the same, the domains are
2920 * obviously the same and it returns 1, at the first difference, it returns 0.
2921 * This is a very fast way to verify this property. It has been shown (with the
2922 * CLooG benchmarks) that operations on equal domains are 17% of all the
2923 * polyhedral computations. For 75% of the actually identical domains, this
2924 * function answer that they are the same and allow to give immediately the
2925 * trivial solution instead of calling the heavy general functions of PolyLib.
2926 * - August 22th 2003: first version.
2927 * - June 21rd 2005: Adaptation for GMP (based on S. Verdoolaege's version of
2931 cloog_domain_lazy_equal (CloogDomain
* d1
, CloogDomain
* d2
)
2934 polyhedra_union u1
= cloog_domain_upol (d1
);
2935 polyhedra_union u2
= cloog_domain_upol (d2
);
2939 if ((cloog_upol_nbc (u1
) != cloog_upol_nbc (u2
)) ||
2940 (cloog_upol_dim (u1
) != cloog_upol_dim (u2
)))
2943 for (i
= 0; i
< (int) cloog_upol_nbc (u1
); i
++)
2944 for (j
= 0; j
< (int) cloog_upol_dim (u1
) + 2; j
++)
2945 if (value_ne (cloog_upol_polyhedron (u1
)->Constraint
[i
][j
],
2946 cloog_upol_polyhedron (u2
)->Constraint
[i
][j
]))
2949 u1
= cloog_upol_next (u1
);
2950 u2
= cloog_upol_next (u2
);
2961 * cloog_domain_lazy_block function:
2962 * This function returns 1 if the two domains d1 and d2 given as input are the
2963 * same (possibly except for a dimension equal to a constant where we accept
2964 * a difference of 1) AND if we are sure that there are no other domain in
2965 * the code generation problem that may put integral points between those of
2966 * d1 and d2 (0 otherwise). In fact this function answers the question "can I
2967 * safely consider the two domains as only one with two statements (a block) ?".
2968 * This function is lazy: it asks for very standard scattering representation
2969 * (only one constraint per dimension which is an equality, and the constraints
2970 * are ordered per dimension depth: the left hand side of the constraint matrix
2971 * is the identity) and will answer NO at the very first problem.
2972 * - d1 and d2 are the two domains to check for blocking,
2973 * - scattering is the linked list of all domains,
2974 * - scattdims is the total number of scattering dimentions.
2976 * - April 30th 2005: beginning
2977 * - June 9th 2005: first working version.
2978 * - June 10th 2005: debugging.
2979 * - June 21rd 2005: Adaptation for GMP.
2980 * - October 16th 2005: (debug) some false blocks have been removed.
2983 cloog_domain_lazy_block (CloogDomain
*d1
, CloogDomain
*d2
, CloogDomainList
*scattering
, int scattdims
)
2985 int i
, j
, difference
= 0, different_constraint
= 0, nbc
;
2987 Value date1
, date2
, date3
, temp
;
2990 /* Some basic checks: we only accept convex domains, with same constraint
2991 * and dimension numbers.
2993 if (!cloog_domain_isconvex (d1
) || !cloog_domain_isconvex (d2
) ||
2994 (cloog_domain_nbconstraints (d1
) != cloog_domain_nbconstraints (d2
)) ||
2995 (cloog_domain_dim (d1
) != cloog_domain_dim (d2
)))
2998 p1
= cloog_domain_polyhedron (d1
);
2999 p2
= cloog_domain_polyhedron (d2
);
3000 nbc
= cloog_domain_nbconstraints (d1
);
3001 dim1
= cloog_domain_dim (d1
);
3002 dim2
= cloog_domain_dim (d2
);
3004 /* There should be only one difference between the two domains, it
3005 * has to be at the constant level and the difference must be of +1,
3006 * moreover, after the difference all domain coefficient has to be 0.
3007 * The matrix shape is:
3009 * |===========|=====|<- 0 line
3010 * |===========|=====|
3011 * |===========|====?|<- different_constraint line (found here)
3012 * |===========|0000=|
3013 * |===========|0000=|<- pX->NbConstraints line
3016 * | | (pX->Dimension + 2) column
3017 * | scattdims column
3021 value_init_c (temp
);
3022 for (i
= 0; i
< nbc
; i
++)
3024 if (difference
== 0)
3025 { /* All elements except scalar must be equal. */
3026 for (j
= 0; j
< dim1
+ 1; j
++)
3027 if (value_ne (p1
->Constraint
[i
][j
],
3028 p2
->Constraint
[i
][j
]))
3030 value_clear_c (temp
);
3033 /* The scalar may differ from +1 (now j=(p1->Dimension + 1)). */
3034 if (value_ne (p1
->Constraint
[i
][j
],
3035 p2
->Constraint
[i
][j
]))
3037 value_increment (temp
, p2
->Constraint
[i
][j
]);
3038 if (value_ne (p1
->Constraint
[i
][j
], temp
))
3040 value_clear_c (temp
);
3046 different_constraint
= i
;
3051 { /* Scattering coefficients must be equal. */
3052 for (j
= 0; j
< (scattdims
+ 1); j
++)
3053 if (value_ne (p1
->Constraint
[i
][j
],
3054 p2
->Constraint
[i
][j
]))
3056 value_clear_c (temp
);
3060 /* Domain coefficients must be 0. */
3061 for (; j
< dim1
+ 1; j
++)
3062 if (value_notzero_p (p1
->Constraint
[i
][j
])
3063 || value_notzero_p (p2
->Constraint
[i
][j
]))
3065 value_clear_c (temp
);
3069 /* Scalar must be equal. */
3070 if (value_ne (p1
->Constraint
[i
][j
],
3071 p2
->Constraint
[i
][j
]))
3073 value_clear_c (temp
);
3078 value_clear_c (temp
);
3080 /* If the domains are exactly the same, this is a block. */
3081 if (difference
== 0)
3084 /* Now a basic check that the constraint with the difference is an
3085 * equality of a dimension with a constant.
3087 for (i
= 0; i
<= different_constraint
; i
++)
3088 if (value_notzero_p (p1
->Constraint
[different_constraint
][i
]))
3091 if (value_notone_p (p1
->Constraint
[different_constraint
][different_constraint
+ 1]))
3094 for (i
= different_constraint
+ 2; i
< dim1
+ 1; i
++)
3095 if (value_notzero_p (p1
->Constraint
[different_constraint
][i
]))
3098 /* For the moment, d1 and d2 are a block candidate. There remains to check
3099 * that there is no other domain that may put an integral point between
3100 * them. In our lazy test we ensure this property by verifying that the
3101 * constraint matrices have a very strict shape: let us consider that the
3102 * dimension with the difference is d. Then the first d dimensions are
3103 * defined in their depth order using equalities (thus the first column begins
3104 * with d zeroes, there is a d*d identity matrix and a zero-matrix for
3105 * the remaining simensions). If a domain can put integral points between the
3106 * domains of the block candidate, this means that the other entries on the
3107 * first d constraints are equal to those of d1 or d2. Thus we are looking for
3108 * such a constraint system, if it exists d1 and d2 is considered to not be
3109 * a block, it is a bock otherwise.
3111 * 1. Only equalities (for the first different_constraint+1 lines).
3112 * | 2. Must be the identity.
3113 * | | 3. Must be zero.
3114 * | | | 4. Elements are equal, the last one is either date1 or 2.
3117 * |0|100|00000|=====|<- 0 line
3118 * |0|010|00000|=====|
3119 * |0|001|00000|====?|<- different_constraint line
3120 * |*|***|*****|*****|
3121 * |*|***|*****|*****|<- pX->NbConstraints line
3124 * | | | (pX->Dimension + 2) column
3125 * | | scattdims column
3126 * | different_constraint+1 column
3130 /* Step 1 and 2. This is only necessary to check one domain because
3131 * we checked that they are equal on this part before.
3133 for (i
= 0; i
<= different_constraint
; i
++)
3135 for (j
= 0; j
< i
+ 1; j
++)
3136 if (value_notzero_p (p1
->Constraint
[i
][j
]))
3139 if (value_notone_p (p1
->Constraint
[i
][i
+ 1]))
3142 for (j
= i
+ 2; j
<= different_constraint
+ 1; j
++)
3143 if (value_notzero_p (p1
->Constraint
[i
][j
]))
3148 for (i
= 0; i
<= different_constraint
; i
++)
3149 for (j
= different_constraint
+ 2; j
<= scattdims
; j
++)
3150 if (value_notzero_p (p1
->Constraint
[i
][j
]))
3153 value_init_c (date1
);
3154 value_init_c (date2
);
3155 value_init_c (date3
);
3157 /* Now we have to check that the two different dates are unique. */
3158 value_assign (date1
, p1
->Constraint
[different_constraint
][dim1
+ 1]);
3159 value_assign (date2
, p2
->Constraint
[different_constraint
][dim2
+ 1]);
3161 /* Step 4. We check all domains except d1 and d2 and we look for at least
3162 * a difference with d1 or d2 on the first different_constraint+1 dimensions.
3164 while (scattering
!= NULL
)
3166 if ((cloog_domain (scattering
) != d1
)
3167 && (cloog_domain (scattering
) != d2
))
3169 CloogDomain
*d3
= cloog_domain (scattering
);
3170 polyhedron p3
= cloog_domain_polyhedron (d3
);
3171 int dim3
= cloog_domain_dim (d3
);
3173 value_assign (date3
,
3174 p3
->Constraint
[different_constraint
][dim3
+ 1]);
3177 if (value_ne (date3
, date2
) && value_ne (date3
, date1
))
3180 for (i
= 0; (i
< different_constraint
) && (!difference
); i
++)
3181 for (j
= 0; (j
< dim3
+ 2) && !difference
; j
++)
3183 (p1
->Constraint
[i
][j
],
3184 p3
->Constraint
[i
][j
]))
3187 for (j
= 0; (j
< dim3
+ 1) && !difference
; j
++)
3189 (p1
->Constraint
[different_constraint
][j
],
3190 p3
->Constraint
[different_constraint
][j
]))
3195 value_clear_c (date1
);
3196 value_clear_c (date2
);
3197 value_clear_c (date3
);
3202 scattering
= cloog_next_domain (scattering
);
3205 value_clear_c (date1
);
3206 value_clear_c (date2
);
3207 value_clear_c (date3
);
3213 * cloog_domain_lazy_disjoint function:
3214 * This function returns 1 if the domains given as input are disjoint, 0 if it
3215 * is unable to decide. This function finds the unknown with fixed values in
3216 * both domains (on a given constraint, their column entry is not zero and
3217 * only the constant coefficient can be different from zero) and verify that
3218 * their values are the same. If not, the domains are obviously disjoint and
3219 * it returns 1, if there is not such case it returns 0. This is a very fast
3220 * way to verify this property. It has been shown (with the CLooG benchmarks)
3221 * that operations on disjoint domains are 36% of all the polyhedral
3222 * computations. For 94% of the actually identical domains, this
3223 * function answer that they are disjoint and allow to give immediately the
3224 * trivial solution instead of calling the heavy general functions of PolyLib.
3225 * - August 22th 2003: first version.
3226 * - June 21rd 2005: Adaptation for GMP (based on S. Verdoolaege's version of
3230 cloog_domain_lazy_disjoint (CloogDomain
* d1
, CloogDomain
* d2
)
3232 int i1
, j1
, i2
, j2
, scat_dim
, nbc1
, nbc2
;
3237 if (!cloog_domain_isconvex (d1
) || !cloog_domain_isconvex (d2
))
3240 p1
= cloog_domain_polyhedron (d1
);
3241 p2
= cloog_domain_polyhedron (d2
);
3242 nbc1
= cloog_domain_nbconstraints (d1
);
3243 nbc2
= cloog_domain_nbconstraints (d2
);
3244 dim1
= cloog_domain_dim (d1
);
3245 dim2
= cloog_domain_dim (d2
);
3246 value_init_c (scat_val
);
3248 for (i1
= 0; i1
< nbc1
; i1
++)
3250 if (value_notzero_p (p1
->Constraint
[i1
][0]))
3254 while (value_zero_p (p1
->Constraint
[i1
][scat_dim
]) &&
3258 if (value_notone_p (p1
->Constraint
[i1
][scat_dim
]))
3262 for (j1
= scat_dim
+ 1; j1
<= dim1
; j1
++)
3263 if (value_notzero_p (p1
->Constraint
[i1
][j1
]))
3269 value_assign (scat_val
,
3270 p1
->Constraint
[i1
][dim1
+ 1]);
3272 for (i2
= 0; i2
< nbc2
; i2
++)
3274 for (j2
= 0; j2
< scat_dim
; j2
++)
3275 if (value_notzero_p (p2
->Constraint
[i2
][j2
]))
3278 if ((j2
!= scat_dim
)
3280 value_notone_p (p2
->Constraint
[i2
][scat_dim
]))
3283 for (j2
= scat_dim
+ 1; j2
< dim2
; j2
++)
3284 if (value_notzero_p (p2
->Constraint
[i2
][j2
]))
3291 (p2
->Constraint
[i2
][dim2
+ 1], scat_val
))
3293 value_clear_c (scat_val
);
3300 value_clear_c (scat_val
);
3306 * cloog_domain_list_lazy_same function:
3307 * This function returns 1 if two domains in the list are the same, 0 if it
3308 * is unable to decide.
3309 * - February 9th 2004: first version.
3312 cloog_domain_list_lazy_same (CloogDomainList
* list
)
3313 { /*int i=1, j=1 ; */
3314 CloogDomainList
*current
, *next
;
3317 while (current
!= NULL
)
3319 next
= cloog_next_domain (current
);
3321 while (next
!= NULL
)
3323 if (cloog_domain_lazy_equal (cloog_domain (current
),
3324 cloog_domain (next
)))
3325 { /*printf("Same domains: %d and %d\n",i,j) ; */
3329 next
= cloog_next_domain (next
);
3332 current
= cloog_next_domain (current
);
3339 * cloog_domain_cut_first function:
3340 * this function returns a CloogDomain structure with everything except the
3341 * first part of the polyhedra union of the input domain as domain. After a call
3342 * to this function, there remains in the CloogDomain structure provided as
3343 * input only the first part of the original polyhedra union.
3344 * - April 20th 2005: first version, extracted from different part of loop.c.
3347 cloog_domain_cut_first (CloogDomain
* domain
)
3351 if (domain
&& cloog_domain_polyhedron (domain
))
3353 if (!cloog_upol_next (cloog_domain_upol (domain
)))
3356 rest
= cloog_new_domain (cloog_upol_next (cloog_domain_upol (domain
)));
3357 cloog_upol_set_next (cloog_domain_upol (domain
), NULL
);
3362 return cloog_check_domain (rest
);
3367 * cloog_domain_lazy_isscalar function:
3368 * this function returns 1 if the dimension 'dimension' in the domain 'domain'
3369 * is scalar, this means that the only constraint on this dimension must have
3370 * the shape "x.dimension + scalar = 0" with x an integral variable. This
3371 * function is lazy since we only accept x=1 (further calculations are easier
3373 * - June 14th 2005: first version.
3374 * - June 21rd 2005: Adaptation for GMP.
3377 cloog_domain_lazy_isscalar (CloogDomain
* domain
, int dimension
)
3380 polyhedron p
= cloog_domain_polyhedron (domain
);
3381 int nbc
= cloog_domain_nbconstraints (domain
);
3382 int dim
= cloog_domain_dim (domain
);
3384 /* For each constraint... */
3385 for (i
= 0; i
< nbc
; i
++)
3386 { /* ...if it is concerned by the potentially scalar dimension... */
3388 (p
->Constraint
[i
][dimension
+ 1]))
3389 { /* ...check that the constraint has the shape "dimension + scalar = 0". */
3390 for (j
= 0; j
<= dimension
; j
++)
3391 if (value_notzero_p (p
->Constraint
[i
][j
]))
3395 (p
->Constraint
[i
][dimension
+ 1]))
3398 for (j
= dimension
+ 2; j
< dim
+ 1; j
++)
3399 if (value_notzero_p (p
->Constraint
[i
][j
]))
3409 * cloog_domain_scalar function:
3410 * when we call this function, we know that "dimension" is a scalar dimension,
3411 * this function finds the scalar value in "domain" and returns it in "value".
3412 * - June 30th 2005: first version.
3415 cloog_domain_scalar (CloogDomain
* domain
, int dimension
, Value
* value
)
3418 polyhedron p
= cloog_domain_polyhedron (domain
);
3419 int nbc
= cloog_domain_nbconstraints (domain
);
3420 int dim
= cloog_domain_dim (domain
);
3422 /* For each constraint... */
3423 for (i
= 0; i
< nbc
; i
++)
3424 { /* ...if it is the equality defining the scalar dimension... */
3426 (p
->Constraint
[i
][dimension
+ 1])
3427 && value_zero_p (p
->Constraint
[i
][0]))
3428 { /* ...Then send the scalar value. */
3429 value_assign (*value
, p
->Constraint
[i
][dim
+ 1]);
3430 value_oppose (*value
, *value
);
3435 /* We should have found a scalar value: if not, there is an error. */
3436 fprintf (stderr
, "[CLooG]ERROR: dimension %d is not scalar as expected.\n",
3443 * cloog_domain_erase_dimension function:
3444 * this function returns a CloogDomain structure builds from 'domain' where
3445 * we removed the dimension 'dimension' and every constraint considering this
3446 * dimension. This is not a projection ! Every data concerning the
3447 * considered dimension is simply erased.
3448 * - June 14th 2005: first version.
3449 * - June 21rd 2005: Adaptation for GMP.
3452 cloog_domain_erase_dimension (CloogDomain
* domain
, int dimension
)
3454 int i
, j
, mi
, nb_dim
, nbc
;
3455 CloogMatrix
*matrix
;
3456 CloogDomain
*erased
;
3457 polyhedron p
= cloog_domain_polyhedron (domain
);
3459 nb_dim
= cloog_domain_dim (domain
);
3460 nbc
= cloog_domain_nbconstraints (domain
);
3462 /* The matrix is one column less and at least one constraint less. */
3463 matrix
= cloog_matrix_alloc (nbc
- 1, nb_dim
+ 1);
3465 /* mi is the constraint counter for the matrix. */
3467 for (i
= 0; i
< nbc
; i
++)
3468 if (value_zero_p (p
->Constraint
[i
][dimension
+ 1]))
3470 for (j
= 0; j
<= dimension
; j
++)
3471 value_assign (matrix
->p
[mi
][j
],
3472 p
->Constraint
[i
][j
]);
3474 for (j
= dimension
+ 2; j
< nb_dim
+ 2; j
++)
3475 value_assign (matrix
->p
[mi
][j
- 1],
3476 p
->Constraint
[i
][j
]);
3481 erased
= cloog_domain_matrix2domain (matrix
);
3482 cloog_matrix_free (matrix
);
3484 return cloog_check_domain (erased
);
3487 /* Number of polyhedra inside the union of disjoint polyhedra. */
3490 cloog_domain_nb_polyhedra (CloogDomain
* domain
)
3493 polyhedra_union upol
= cloog_domain_upol (domain
);
3498 upol
= cloog_upol_next (upol
);
3505 cloog_domain_print_polyhedra (FILE * foo
, CloogDomain
* domain
)
3507 polyhedra_union upol
= cloog_domain_upol (domain
);
3509 while (upol
!= NULL
)
3511 CloogMatrix
*matrix
= cloog_upol_domain2matrix (upol
);
3512 cloog_matrix_print (foo
, matrix
);
3513 cloog_matrix_free (matrix
);
3514 upol
= cloog_upol_next (upol
);
3519 debug_cloog_domain (CloogDomain
*domain
)
3521 cloog_domain_print_polyhedra (stderr
, domain
);
3525 debug_cloog_matrix (CloogMatrix
*m
)
3527 cloog_matrix_print (stderr
, m
);
3531 debug_value (Value v
)
3533 value_print (stderr
, VALUE_FMT
, v
);
3537 debug_values (Value
*p
, int length
)
3540 for (i
= 0; i
< length
; i
++)
3544 polyhedra_union
cloog_new_upol (polyhedron p
)
3546 polyhedra_union ppl
=
3547 (polyhedra_union
) malloc (sizeof (struct polyhedra_union
));
3548 ppl
->_polyhedron
= p
;
3553 Vector
*Vector_Alloc (unsigned length
)
3556 Vector
*vector
= (Vector
*) malloc (sizeof (Vector
));
3558 vector
->Size
= length
;
3559 vector
->p
= (Value
*) malloc (length
* sizeof (Value
));
3561 for (i
= 0; i
< length
; i
++)
3562 value_init (vector
->p
[i
]);
3567 polyhedron
cloog_new_pol (int dim
, int nrows
)
3570 polyhedron res
= (polyhedron
) malloc (sizeof (struct polyhedron
));
3571 int ncolumns
= dim
+ 2;
3572 int n
= nrows
* ncolumns
;
3573 Value
*p
= (Value
*) malloc (n
* sizeof (Value
));
3575 res
->Dimension
= dim
;
3576 res
->NbConstraints
= nrows
;
3577 res
->Constraint
= (Value
**) malloc (nrows
* sizeof (Value
*));
3579 for (i
= 0; i
< n
; ++i
)
3582 for (i
= 0; i
< nrows
; i
++, p
+= ncolumns
)
3583 res
->Constraint
[i
] = p
;