isl_obj_list_print: properly return isl_printer
[barvinok.git] / evalue_isl.c
blob65c29f7fad5117f13178de3af14fba9281816b7c
1 #include <isl_set_polylib.h>
2 #include <isl/constraint.h>
3 #include <isl/seq.h>
4 #include <barvinok/evalue.h>
6 static __isl_give isl_qpolynomial *extract_base(__isl_take isl_dim *dim,
7 const evalue *e)
9 int i;
10 isl_ctx *ctx;
11 isl_vec *v;
12 isl_div *div;
13 isl_qpolynomial *base, *c;
14 unsigned nparam;
16 if (!dim)
17 return NULL;
19 if (e->x.p->type == polynomial)
20 return isl_qpolynomial_var(dim, isl_dim_param, e->x.p->pos - 1);
22 ctx = isl_dim_get_ctx(dim);
23 nparam = isl_dim_size(dim, isl_dim_param);
24 v = isl_vec_alloc(ctx, 2 + nparam);
25 if (!v)
26 goto error;
28 isl_seq_clr(v->el + 2, nparam);
29 evalue_extract_affine(&e->x.p->arr[0], v->el + 2, &v->el[1], &v->el[0]);
31 div = isl_div_alloc(isl_dim_copy(dim));
32 isl_div_set_constant(div, v->el[1]);
33 isl_div_set_denominator(div, v->el[0]);
35 for (i = 0; i < nparam; ++i)
36 isl_div_set_coefficient(div, isl_dim_param, i, v->el[2 + i]);
38 base = isl_qpolynomial_div(div);
40 if (e->x.p->type == fractional) {
41 base = isl_qpolynomial_neg(base);
43 c = isl_qpolynomial_rat_cst(isl_dim_copy(dim), v->el[1], v->el[0]);
44 base = isl_qpolynomial_add(base, c);
46 for (i = 0; i < nparam; ++i) {
47 isl_qpolynomial *t;
48 c = isl_qpolynomial_rat_cst(isl_dim_copy(dim),
49 v->el[2 + i], v->el[0]);
50 t = isl_qpolynomial_var(isl_dim_copy(dim),
51 isl_dim_param, i);
52 t = isl_qpolynomial_mul(c, t);
53 base = isl_qpolynomial_add(base, t);
57 isl_vec_free(v);
58 isl_dim_free(dim);
60 return base;
61 error:
62 isl_dim_free(dim);
63 return NULL;
66 static int type_offset(enode *p)
68 return p->type == fractional ? 1 :
69 p->type == flooring ? 1 : 0;
72 __isl_give isl_qpolynomial *isl_qpolynomial_from_evalue(__isl_take isl_dim *dim,
73 const evalue *e)
75 int i;
76 isl_qpolynomial *qp;
77 isl_qpolynomial *base;
78 int offset;
80 if (EVALUE_IS_NAN(*e))
81 return isl_qpolynomial_infty(dim);
82 if (value_notzero_p(e->d))
83 return isl_qpolynomial_rat_cst(dim, e->x.n, e->d);
85 offset = type_offset(e->x.p);
87 assert(e->x.p->type == polynomial ||
88 e->x.p->type == flooring ||
89 e->x.p->type == fractional);
90 assert(e->x.p->size >= 1 + offset);
92 base = extract_base(isl_dim_copy(dim), e);
93 qp = isl_qpolynomial_from_evalue(isl_dim_copy(dim),
94 &e->x.p->arr[e->x.p->size - 1]);
96 for (i = e->x.p->size - 2; i >= offset; --i) {
97 qp = isl_qpolynomial_mul(qp, isl_qpolynomial_copy(base));
98 qp = isl_qpolynomial_add(qp,
99 isl_qpolynomial_from_evalue(isl_dim_copy(dim),
100 &e->x.p->arr[i]));
103 isl_qpolynomial_free(base);
104 isl_dim_free(dim);
106 return qp;
109 static __isl_give isl_pw_qpolynomial *guarded_evalue2pwqp(__isl_take isl_set *set,
110 const evalue *e);
112 static __isl_give isl_pw_qpolynomial *relation2pwqp(__isl_take isl_set *set,
113 const evalue *e)
115 int i;
116 isl_vec *vec;
117 isl_dim *dim;
118 isl_ctx *ctx;
119 unsigned nparam;
120 isl_pw_qpolynomial *pwqp;
121 struct isl_constraint *c;
122 struct isl_basic_set *bset;
123 struct isl_set *guard;
124 const evalue *fract;
126 if (!set || !e)
127 goto error;
129 if (e->x.p->size == 1) {
130 dim = isl_set_get_dim(set);
131 isl_set_free(set);
132 return isl_pw_qpolynomial_zero(dim);
135 ctx = isl_set_get_ctx(set);
136 isl_assert(ctx, e->x.p->size > 0, goto error);
137 isl_assert(ctx, e->x.p->size <= 3, goto error);
138 isl_assert(ctx, value_zero_p(e->x.p->arr[0].d), goto error);
139 isl_assert(ctx, e->x.p->arr[0].x.p->type == fractional, goto error);
140 fract = &e->x.p->arr[0];
142 nparam = isl_set_dim(set, isl_dim_param);
143 vec = isl_vec_alloc(ctx, 2 + nparam + 1);
144 if (!vec)
145 goto error;
147 isl_seq_clr(vec->el + 2, nparam);
148 evalue_extract_affine(&fract->x.p->arr[0],
149 vec->el + 2, &vec->el[1], &vec->el[0]);
151 dim = isl_set_get_dim(set);
152 dim = isl_dim_add(dim, isl_dim_set, 1);
154 bset = isl_basic_set_universe(dim);
155 c = isl_equality_alloc(isl_dim_copy(dim));
156 isl_int_neg(vec->el[0], vec->el[0]);
157 isl_constraint_set_coefficient(c, isl_dim_set, 0, vec->el[0]);
158 isl_constraint_set_constant(c, vec->el[1]);
159 for (i = 0; i < nparam; ++i)
160 isl_constraint_set_coefficient(c, isl_dim_param, i, vec->el[2+i]);
161 bset = isl_basic_set_add_constraint(bset, c);
162 bset = isl_basic_set_project_out(bset, isl_dim_set, 0, 1);
163 guard = isl_set_from_basic_set(bset);
164 isl_vec_free(vec);
166 pwqp = guarded_evalue2pwqp(isl_set_intersect(isl_set_copy(set),
167 isl_set_copy(guard)),
168 &e->x.p->arr[1]);
170 if (e->x.p->size == 3) {
171 isl_pw_qpolynomial *pwqpc;
172 guard = isl_set_complement(guard);
173 pwqpc = guarded_evalue2pwqp(isl_set_intersect(isl_set_copy(set),
174 isl_set_copy(guard)),
175 &e->x.p->arr[2]);
176 pwqp = isl_pw_qpolynomial_add_disjoint(pwqp, pwqpc);
179 isl_set_free(set);
180 isl_set_free(guard);
182 return pwqp;
183 error:
184 isl_set_free(set);
185 return NULL;
188 static __isl_give isl_pw_qpolynomial *guarded_evalue2pwqp(__isl_take isl_set *set,
189 const evalue *e)
191 isl_qpolynomial *qp;
193 if (value_zero_p(e->d) && e->x.p->type == relation)
194 return relation2pwqp(set, e);
196 qp = isl_qpolynomial_from_evalue(isl_set_get_dim(set), e);
198 return isl_pw_qpolynomial_alloc(set, qp);
201 __isl_give isl_pw_qpolynomial *isl_pw_qpolynomial_from_evalue(__isl_take isl_dim *dim, const evalue *e)
203 int i;
204 isl_pw_qpolynomial *pwqp;
206 if (!dim || !e)
207 goto error;
208 if (EVALUE_IS_ZERO(*e))
209 return isl_pw_qpolynomial_zero(dim);
211 if (value_notzero_p(e->d)) {
212 isl_set *set = isl_set_universe(isl_dim_copy(dim));
213 isl_qpolynomial *qp = isl_qpolynomial_rat_cst(dim, e->x.n, e->d);
214 return isl_pw_qpolynomial_alloc(set, qp);
217 assert(!EVALUE_IS_NAN(*e));
219 assert(e->x.p->type == partition);
221 pwqp = isl_pw_qpolynomial_zero(isl_dim_copy(dim));
223 for (i = 0; i < e->x.p->size/2; ++i) {
224 Polyhedron *D = EVALUE_DOMAIN(e->x.p->arr[2 * i]);
225 isl_set *set = isl_set_new_from_polylib(D, isl_dim_copy(dim));
226 isl_pw_qpolynomial *pwqp_i;
228 pwqp_i = guarded_evalue2pwqp(set, &e->x.p->arr[2 * i + 1]);
230 pwqp = isl_pw_qpolynomial_add_disjoint(pwqp, pwqp_i);
233 isl_dim_free(dim);
235 return pwqp;
236 error:
237 isl_dim_free(dim);
238 return NULL;
241 static evalue *evalue_pow(evalue *e, int exp)
243 evalue *pow;
245 if (exp == 1)
246 return e;
248 pow = evalue_dup(e);
249 while (--exp > 0)
250 emul(e, pow);
252 evalue_free(e);
254 return pow;
257 static evalue *div2evalue(__isl_take isl_div *div)
259 int i;
260 isl_ctx *ctx;
261 isl_vec *vec;
262 unsigned dim;
263 unsigned nparam;
264 evalue *e;
265 evalue *aff;
267 if (!div)
268 return NULL;
270 dim = isl_div_dim(div, isl_dim_set);
271 nparam = isl_div_dim(div, isl_dim_param);
273 ctx = isl_div_get_ctx(div);
274 vec = isl_vec_alloc(ctx, 1 + dim + nparam + 1);
275 if (!vec)
276 goto error;
277 for (i = 0; i < dim; ++i)
278 isl_div_get_coefficient(div, isl_dim_set, i, &vec->el[1 + i]);
279 for (i = 0; i < nparam; ++i)
280 isl_div_get_coefficient(div, isl_dim_param, i,
281 &vec->el[1 + dim + i]);
282 isl_div_get_denominator(div, &vec->el[0]);
283 isl_div_get_constant(div, &vec->el[1 + dim + nparam]);
285 e = isl_alloc_type(ctx, evalue);
286 if (!e)
287 goto error;
288 value_init(e->d);
289 value_set_si(e->d, 0);
290 e->x.p = new_enode(flooring, 3, -1);
291 evalue_set_si(&e->x.p->arr[1], 0, 1);
292 evalue_set_si(&e->x.p->arr[2], 1, 1);
293 value_clear(e->x.p->arr[0].d);
294 aff = affine2evalue(vec->el + 1, vec->el[0], dim + nparam);
295 e->x.p->arr[0] = *aff;
296 free(aff);
297 isl_vec_free(vec);
298 isl_div_free(div);
299 return e;
300 error:
301 isl_vec_free(vec);
302 isl_div_free(div);
303 return NULL;
306 static int add_term(__isl_take isl_term *term, void *user)
308 int i;
309 evalue *sum = (evalue *)user;
310 unsigned nparam;
311 unsigned dim;
312 unsigned n_div;
313 isl_ctx *ctx;
314 isl_int n, d;
315 evalue *e;
317 if (!term)
318 return -1;
320 nparam = isl_term_dim(term, isl_dim_param);
321 dim = isl_term_dim(term, isl_dim_set);
322 n_div = isl_term_dim(term, isl_dim_div);
324 ctx = isl_term_get_ctx(term);
325 e = isl_alloc_type(ctx, evalue);
326 if (!e)
327 goto error;
329 isl_int_init(n);
330 isl_int_init(d);
332 isl_term_get_num(term, &n);
333 isl_term_get_den(term, &d);
334 value_init(e->d);
335 evalue_set(e, n, d);
337 for (i = 0; i < dim; ++i) {
338 evalue *pow;
339 int exp = isl_term_get_exp(term, isl_dim_set, i);
341 if (!exp)
342 continue;
344 pow = evalue_pow(evalue_var(i), exp);
345 emul(pow, e);
346 evalue_free(pow);
349 for (i = 0; i < nparam; ++i) {
350 evalue *pow;
351 int exp = isl_term_get_exp(term, isl_dim_param, i);
353 if (!exp)
354 continue;
356 pow = evalue_pow(evalue_var(dim + i), exp);
357 emul(pow, e);
358 evalue_free(pow);
361 for (i = 0; i < n_div; ++i) {
362 evalue *pow;
363 evalue *floor;
364 isl_div *div;
365 int exp = isl_term_get_exp(term, isl_dim_div, i);
367 if (!exp)
368 continue;
370 div = isl_term_get_div(term, i);
371 floor = div2evalue(div);
372 pow = evalue_pow(floor, exp);
373 emul(pow, e);
374 evalue_free(pow);
377 eadd(e, sum);
378 evalue_free(e);
380 isl_int_clear(n);
381 isl_int_clear(d);
383 isl_term_free(term);
385 return 0;
386 error:
387 isl_term_free(term);
388 return -1;
391 evalue *isl_qpolynomial_to_evalue(__isl_keep isl_qpolynomial *qp)
393 evalue *e;
395 e = evalue_zero();
396 if (!e)
397 return NULL;
399 if (isl_qpolynomial_foreach_term(qp, add_term, e) < 0)
400 goto error;
402 return e;
403 error:
404 evalue_free(e);
405 return NULL;
408 static int add_guarded_qp(__isl_take isl_set *set, __isl_take isl_qpolynomial *qp,
409 void *user)
411 Polyhedron *D;
412 evalue *e = NULL;
413 evalue *f;
414 evalue *sum = (evalue *)user;
415 unsigned dim;
417 e = isl_alloc_type(isl_set_get_ctx(set), evalue);
418 if (!e)
419 goto error;
421 D = isl_set_to_polylib(set);
422 if (!D)
423 goto error;
425 f = isl_qpolynomial_to_evalue(qp);
426 if (!e) {
427 Domain_Free(D);
428 goto error;
431 dim = isl_set_dim(set, isl_dim_param) + isl_set_dim(set, isl_dim_set);
432 value_init(e->d);
433 e->x.p = new_enode(partition, 2, D->Dimension);
434 EVALUE_SET_DOMAIN(e->x.p->arr[0], D);
436 value_clear(e->x.p->arr[1].d);
437 e->x.p->arr[1] = *f;
438 free(f);
440 eadd(e, sum);
441 evalue_free(e);
443 isl_set_free(set);
444 isl_qpolynomial_free(qp);
446 return 0;
447 error:
448 free(e);
449 isl_set_free(set);
450 isl_qpolynomial_free(qp);
451 return -1;
454 evalue *isl_pw_qpolynomial_to_evalue(__isl_keep isl_pw_qpolynomial *pwqp)
456 evalue *e;
458 if (!pwqp)
459 return NULL;
460 e = evalue_zero();
462 if (isl_pw_qpolynomial_foreach_piece(pwqp, add_guarded_qp, e))
463 goto error;
465 return e;
466 error:
467 evalue_free(e);
468 return NULL;