4 #include "isl_piplib.h"
5 #include "isl_sample_piplib.h"
7 static void swap_inequality(struct isl_basic_set
*bset
, int a
, int b
)
9 isl_int
*t
= bset
->ineq
[a
];
10 bset
->ineq
[a
] = bset
->ineq
[b
];
14 static struct isl_mat
*independent_bounds(struct isl_ctx
*ctx
,
15 struct isl_basic_set
*bset
)
18 struct isl_mat
*dirs
= NULL
;
19 struct isl_mat
*bounds
= NULL
;
24 bounds
= isl_mat_alloc(ctx
, 1+bset
->dim
, 1+bset
->dim
);
28 isl_int_set_si(bounds
->row
[0][0], 1);
29 isl_seq_clr(bounds
->row
[0]+1, bset
->dim
);
32 if (bset
->n_ineq
== 0)
35 dirs
= isl_mat_alloc(ctx
, bset
->dim
, bset
->dim
);
37 isl_mat_free(ctx
, bounds
);
40 isl_seq_cpy(dirs
->row
[0], bset
->ineq
[0]+1, dirs
->n_col
);
41 isl_seq_cpy(bounds
->row
[1], bset
->ineq
[0], bounds
->n_col
);
42 for (j
= 1, n
= 1; n
< bset
->dim
&& j
< bset
->n_ineq
; ++j
) {
45 isl_seq_cpy(dirs
->row
[n
], bset
->ineq
[j
]+1, dirs
->n_col
);
47 pos
= isl_seq_first_non_zero(dirs
->row
[n
], dirs
->n_col
);
50 for (i
= 0; i
< n
; ++i
) {
52 pos_i
= isl_seq_first_non_zero(dirs
->row
[i
], dirs
->n_col
);
57 isl_seq_elim(dirs
->row
[n
], dirs
->row
[i
], pos
,
59 pos
= isl_seq_first_non_zero(dirs
->row
[n
], dirs
->n_col
);
67 isl_int
*t
= dirs
->row
[n
];
68 for (k
= n
; k
> i
; --k
)
69 dirs
->row
[k
] = dirs
->row
[k
-1];
73 isl_seq_cpy(bounds
->row
[n
], bset
->ineq
[j
], bounds
->n_col
);
75 isl_mat_free(ctx
, dirs
);
80 /* Skew into positive orthant and project out lineality space */
81 static struct isl_basic_set
*isl_basic_set_skew_to_positive_orthant(
82 struct isl_basic_set
*bset
, struct isl_mat
**T
)
84 struct isl_mat
*U
= NULL
;
85 struct isl_mat
*bounds
= NULL
;
87 unsigned old_dim
, new_dim
;
95 isl_assert(ctx
, bset
->nparam
== 0, goto error
);
96 isl_assert(ctx
, bset
->n_div
== 0, goto error
);
97 isl_assert(ctx
, bset
->n_eq
== 0, goto error
);
99 /* Try to move (multiples of) unit rows up. */
100 for (i
= 0, j
= 0; i
< bset
->n_ineq
; ++i
) {
101 int pos
= isl_seq_first_non_zero(bset
->ineq
[i
]+1,
105 if (isl_seq_first_non_zero(bset
->ineq
[i
]+1+pos
+1,
106 bset
->dim
-pos
-1) >= 0)
109 swap_inequality(bset
, i
, j
);
112 bounds
= independent_bounds(ctx
, bset
);
116 new_dim
= bounds
->n_row
- 1;
117 bounds
= isl_mat_left_hermite(ctx
, bounds
, 1, &U
, NULL
);
120 U
= isl_mat_drop_cols(ctx
, U
, 1 + new_dim
, old_dim
- new_dim
);
121 bset
= isl_basic_set_preimage(ctx
, bset
, isl_mat_copy(ctx
, U
));
125 isl_mat_free(ctx
, bounds
);
128 isl_mat_free(ctx
, bounds
);
129 isl_mat_free(ctx
, U
);
130 isl_basic_set_free(bset
);
134 struct isl_vec
*isl_pip_basic_set_sample(struct isl_basic_set
*bset
)
136 PipOptions
*options
= NULL
;
137 PipMatrix
*domain
= NULL
;
138 PipQuast
*sol
= NULL
;
139 struct isl_vec
*vec
= NULL
;
147 isl_assert(ctx
, bset
->nparam
== 0, goto error
);
148 isl_assert(ctx
, bset
->n_div
== 0, goto error
);
149 bset
= isl_basic_set_skew_to_positive_orthant(bset
, &T
);
153 domain
= isl_basic_map_to_pip((struct isl_basic_map
*)bset
, 0, 0, 0);
157 options
= pip_options_init();
160 sol
= pip_solve(domain
, NULL
, -1, options
);
164 vec
= isl_vec_alloc(ctx
, 0);
165 isl_mat_free(ctx
, T
);
169 vec
= isl_vec_alloc(ctx
, 1 + dim
);
172 isl_int_set_si(vec
->block
.data
[0], 1);
173 for (i
= 0, l
= sol
->list
; l
&& i
< dim
; ++i
, l
= l
->next
) {
174 isl_seq_cpy_from_pip(&vec
->block
.data
[1+i
],
175 &l
->vector
->the_vector
[0], 1);
176 isl_assert(ctx
, !entier_zero_p(l
->vector
->the_deno
[0]),
179 isl_assert(ctx
, i
== dim
, goto error
);
180 vec
= isl_mat_vec_product(ctx
, T
, vec
);
184 pip_options_free(options
);
185 pip_matrix_free(domain
);
187 isl_basic_set_free(bset
);
190 isl_vec_free(ctx
, vec
);
191 isl_basic_set_free(bset
);
195 pip_matrix_free(domain
);