5 #include <barvinok/genfun.h>
6 #include <barvinok/barvinok.h>
7 #include "conversion.h"
8 #include "genfun_constructor.h"
17 bool short_rat_lex_smaller_denominator::operator()(const short_rat
* r1
,
18 const short_rat
* r2
) const
20 return lex_cmp(r1
->d
.power
, r2
->d
.power
) < 0;
23 static void lex_order_terms(struct short_rat
* rat
)
25 for (int i
= 0; i
< rat
->n
.power
.NumRows(); ++i
) {
27 for (int j
= i
+1; j
< rat
->n
.power
.NumRows(); ++j
)
28 if (lex_cmp(rat
->n
.power
[j
], rat
->n
.power
[m
]) < 0)
31 vec_ZZ tmp
= rat
->n
.power
[m
];
32 rat
->n
.power
[m
] = rat
->n
.power
[i
];
33 rat
->n
.power
[i
] = tmp
;
34 QQ tmp_coeff
= rat
->n
.coeff
[m
];
35 rat
->n
.coeff
[m
] = rat
->n
.coeff
[i
];
36 rat
->n
.coeff
[i
] = tmp_coeff
;
41 short_rat::short_rat(Value c
)
44 value2zz(c
, n
.coeff
[0].n
);
46 n
.power
.SetDims(1, 0);
47 d
.power
.SetDims(0, 0);
50 short_rat::short_rat(const QQ
& c
, const vec_ZZ
& num
, const mat_ZZ
& den
)
56 n
.power
.SetDims(1, num
.length());
62 void short_rat::normalize()
64 /* Make all powers in denominator lexico-positive */
65 for (int i
= 0; i
< d
.power
.NumRows(); ++i
) {
67 for (j
= 0; j
< d
.power
.NumCols(); ++j
)
68 if (d
.power
[i
][j
] != 0)
70 assert(j
< d
.power
.NumCols());
71 if (d
.power
[i
][j
] < 0) {
72 d
.power
[i
] = -d
.power
[i
];
73 for (int k
= 0; k
< n
.coeff
.length(); ++k
) {
74 n
.coeff
[k
].n
= -n
.coeff
[k
].n
;
75 n
.power
[k
] += d
.power
[i
];
80 /* Order powers in denominator */
81 lex_order_rows(d
.power
);
84 void short_rat::add(short_rat
*r
)
86 for (int i
= 0; i
< r
->n
.power
.NumRows(); ++i
) {
87 int len
= n
.coeff
.length();
89 for (j
= 0; j
< len
; ++j
)
90 if (r
->n
.power
[i
] == n
.power
[j
])
93 n
.coeff
[j
] += r
->n
.coeff
[i
];
94 if (n
.coeff
[j
].n
== 0) {
96 n
.power
[j
] = n
.power
[len
-1];
97 n
.coeff
[j
] = n
.coeff
[len
-1];
99 int dim
= n
.power
.NumCols();
100 n
.coeff
.SetLength(len
-1);
101 n
.power
.SetDims(len
-1, dim
);
104 int dim
= n
.power
.NumCols();
105 n
.coeff
.SetLength(len
+1);
106 n
.power
.SetDims(len
+1, dim
);
107 n
.coeff
[len
] = r
->n
.coeff
[i
];
108 n
.power
[len
] = r
->n
.power
[i
];
113 bool short_rat::reduced()
115 int dim
= n
.power
.NumCols();
116 lex_order_terms(this);
117 if (n
.power
.NumRows() % 2 == 0) {
118 if (n
.coeff
[0].n
== -n
.coeff
[1].n
&&
119 n
.coeff
[0].d
== n
.coeff
[1].d
) {
120 vec_ZZ step
= n
.power
[1] - n
.power
[0];
122 for (k
= 1; k
< n
.power
.NumRows()/2; ++k
) {
123 if (n
.coeff
[2*k
].n
!= -n
.coeff
[2*k
+1].n
||
124 n
.coeff
[2*k
].d
!= n
.coeff
[2*k
+1].d
)
126 if (step
!= n
.power
[2*k
+1] - n
.power
[2*k
])
129 if (k
== n
.power
.NumRows()/2) {
130 for (k
= 0; k
< d
.power
.NumRows(); ++k
)
131 if (d
.power
[k
] == step
)
133 if (k
< d
.power
.NumRows()) {
134 for (++k
; k
< d
.power
.NumRows(); ++k
)
135 d
.power
[k
-1] = d
.power
[k
];
136 d
.power
.SetDims(k
-1, dim
);
137 for (k
= 1; k
< n
.power
.NumRows()/2; ++k
) {
138 n
.coeff
[k
] = n
.coeff
[2*k
];
139 n
.power
[k
] = n
.power
[2*k
];
141 n
.coeff
.SetLength(k
);
142 n
.power
.SetDims(k
, dim
);
151 gen_fun::gen_fun(Value c
)
153 short_rat
*r
= new short_rat(c
);
154 context
= Universe_Polyhedron(0);
158 void gen_fun::add(const QQ
& c
, const vec_ZZ
& num
, const mat_ZZ
& den
)
163 short_rat
* r
= new short_rat(c
, num
, den
);
165 short_rat_list::iterator i
= term
.find(r
);
166 while (i
!= term
.end()) {
168 if ((*i
)->n
.coeff
.length() == 0) {
171 } else if ((*i
)->reduced()) {
173 /* we've modified term[i], so remove it
174 * and add it back again
188 void gen_fun::add(const QQ
& c
, const gen_fun
*gf
)
191 for (short_rat_list::iterator i
= gf
->term
.begin(); i
!= gf
->term
.end(); ++i
) {
192 for (int j
= 0; j
< (*i
)->n
.power
.NumRows(); ++j
) {
194 p
*= (*i
)->n
.coeff
[j
];
195 add(p
, (*i
)->n
.power
[j
], (*i
)->d
.power
);
200 static void split_param_compression(Matrix
*CP
, mat_ZZ
& map
, vec_ZZ
& offset
)
202 Matrix
*T
= Transpose(CP
);
203 matrix2zz(T
, map
, T
->NbRows
-1, T
->NbColumns
-1);
204 values2zz(T
->p
[T
->NbRows
-1], offset
, T
->NbColumns
-1);
209 * Perform the substitution specified by CP
211 * CP is a homogeneous matrix that maps a set of "compressed parameters"
212 * to the original set of parameters.
214 * This function is applied to a gen_fun computed with the compressed parameters
215 * and adapts it to refer to the original parameters.
217 * That is, if y are the compressed parameters and x = A y + b are the original
218 * parameters, then we want the coefficient of the monomial t^y in the original
219 * generating function to be the coefficient of the monomial u^x in the resulting
220 * generating function.
221 * The original generating function has the form
223 * a t^m/(1-t^n) = a t^m + a t^{m+n} + a t^{m+2n} + ...
225 * Since each term t^y should correspond to a term u^x, with x = A y + b, we want
227 * a u^{A m + b} + a u^{A (m+n) + b} + a u^{A (m+2n) +b} + ... =
229 * = a u^{A m + b}/(1-u^{A n})
231 * Therefore, we multiply the powers m and n in both numerator and denominator by A
232 * and add b to the power in the numerator.
233 * Since the above powers are stored as row vectors m^T and n^T,
234 * we compute, say, m'^T = m^T A^T to obtain m' = A m.
236 * The pair (map, offset) contains the same information as CP.
237 * map is the transpose of the linear part of CP, while offset is the constant part.
239 void gen_fun::substitute(Matrix
*CP
)
243 split_param_compression(CP
, map
, offset
);
244 Polyhedron
*C
= Polyhedron_Image(context
, CP
, 0);
245 Polyhedron_Free(context
);
247 for (short_rat_list::iterator i
= term
.begin(); i
!= term
.end(); ++i
) {
248 (*i
)->d
.power
*= map
;
249 (*i
)->n
.power
*= map
;
250 for (int j
= 0; j
< (*i
)->n
.power
.NumRows(); ++j
)
251 (*i
)->n
.power
[j
] += offset
;
257 vector
<pair
<Vector
*, QQ
> > vertices
;
258 cone(int *pos
) : pos(pos
) {}
261 #ifndef HAVE_COMPRESS_PARMS
262 static Matrix
*compress_parms(Matrix
*M
, unsigned nparam
)
268 struct parallel_polytopes
{
276 barvinok_options
*options
;
278 parallel_polytopes(int n
, Polyhedron
*context
, int nparam
,
279 barvinok_options
*options
) :
280 context(context
), dim(-1), nparam(nparam
),
287 bool add(const QQ
& c
, Polyhedron
*P
) {
290 for (i
= 0; i
< P
->NbEq
; ++i
)
291 if (First_Non_Zero(P
->Constraint
[i
]+1,
292 P
->Dimension
-nparam
) == -1)
297 Polyhedron
*Q
= remove_equalities_p(Polyhedron_Copy(P
), P
->Dimension
-nparam
,
299 POL_ENSURE_VERTICES(Q
);
309 M
= Matrix_Alloc(Q
->NbEq
, Q
->Dimension
+2);
310 Vector_Copy(Q
->Constraint
[0], M
->p
[0], Q
->NbEq
* (Q
->Dimension
+2));
311 CP
= compress_parms(M
, nparam
);
312 T
= align_matrix(CP
, Q
->Dimension
+1);
315 R
= Polyhedron_Preimage(Q
, T
, options
->MaxRays
);
317 Q
= remove_equalities_p(R
, R
->Dimension
-nparam
, NULL
);
319 assert(Q
->NbEq
== 0);
321 if (First_Non_Zero(Q
->Constraint
[Q
->NbConstraints
-1]+1, Q
->Dimension
) == -1)
326 red
= gf_base::create(Polyhedron_Copy(context
), dim
, nparam
, options
);
328 Constraints
= Matrix_Alloc(Q
->NbConstraints
, Q
->Dimension
);
329 for (int i
= 0; i
< Q
->NbConstraints
; ++i
) {
330 Vector_Copy(Q
->Constraint
[i
]+1, Constraints
->p
[i
], Q
->Dimension
);
333 assert(Q
->Dimension
== dim
);
334 for (int i
= 0; i
< Q
->NbConstraints
; ++i
) {
336 for (j
= 0; j
< Constraints
->NbRows
; ++j
)
337 if (Vector_Equal(Q
->Constraint
[i
]+1, Constraints
->p
[j
],
340 assert(j
< Constraints
->NbRows
);
344 for (int i
= 0; i
< Q
->NbRays
; ++i
) {
345 if (!value_pos_p(Q
->Ray
[i
][dim
+1]))
348 Polyhedron
*C
= supporting_cone(Q
, i
);
350 if (First_Non_Zero(C
->Constraint
[C
->NbConstraints
-1]+1,
354 int *pos
= new int[1+C
->NbConstraints
];
355 pos
[0] = C
->NbConstraints
;
357 for (int k
= 0; k
< Constraints
->NbRows
; ++k
) {
358 for (int j
= 0; j
< C
->NbConstraints
; ++j
) {
359 if (Vector_Equal(C
->Constraint
[j
]+1, Constraints
->p
[k
],
366 assert(l
== C
->NbConstraints
);
369 for (j
= 0; j
< cones
.size(); ++j
)
370 if (!memcmp(pos
, cones
[j
].pos
, (1+C
->NbConstraints
)*sizeof(int)))
372 if (j
== cones
.size())
373 cones
.push_back(cone(pos
));
380 for (k
= 0; k
< cones
[j
].vertices
.size(); ++k
)
381 if (Vector_Equal(Q
->Ray
[i
]+1, cones
[j
].vertices
[k
].first
->p
,
385 if (k
== cones
[j
].vertices
.size()) {
386 Vector
*vertex
= Vector_Alloc(Q
->Dimension
+1);
387 Vector_Copy(Q
->Ray
[i
]+1, vertex
->p
, Q
->Dimension
+1);
388 cones
[j
].vertices
.push_back(pair
<Vector
*,QQ
>(vertex
, c
));
390 cones
[j
].vertices
[k
].second
+= c
;
391 if (cones
[j
].vertices
[k
].second
.n
== 0) {
392 int size
= cones
[j
].vertices
.size();
393 Vector_Free(cones
[j
].vertices
[k
].first
);
395 cones
[j
].vertices
[k
] = cones
[j
].vertices
[size
-1];
396 cones
[j
].vertices
.pop_back();
407 for (int i
= 0; i
< cones
.size(); ++i
) {
408 Matrix
*M
= Matrix_Alloc(cones
[i
].pos
[0], 1+Constraints
->NbColumns
+1);
410 for (int j
= 0; j
<cones
[i
].pos
[0]; ++j
) {
411 value_set_si(M
->p
[j
][0], 1);
412 Vector_Copy(Constraints
->p
[cones
[i
].pos
[1+j
]], M
->p
[j
]+1,
413 Constraints
->NbColumns
);
415 Cone
= Constraints2Polyhedron(M
, options
->MaxRays
);
417 for (int j
= 0; j
< cones
[i
].vertices
.size(); ++j
) {
418 red
->base
->do_vertex_cone(cones
[i
].vertices
[j
].second
,
419 Polyhedron_Copy(Cone
),
420 cones
[i
].vertices
[j
].first
->p
, options
);
422 Polyhedron_Free(Cone
);
425 red
->gf
->substitute(CP
);
428 void print(std::ostream
& os
) const {
429 for (int i
= 0; i
< cones
.size(); ++i
) {
431 for (int j
= 0; j
< cones
[i
].pos
[0]; ++j
) {
434 os
<< cones
[i
].pos
[1+j
];
437 for (int j
= 0; j
< cones
[i
].vertices
.size(); ++j
) {
438 Vector_Print(stderr
, P_VALUE_FMT
, cones
[i
].vertices
[j
].first
);
439 os
<< cones
[i
].vertices
[j
].second
<< endl
;
443 ~parallel_polytopes() {
444 for (int i
= 0; i
< cones
.size(); ++i
) {
445 delete [] cones
[i
].pos
;
446 for (int j
= 0; j
< cones
[i
].vertices
.size(); ++j
)
447 Vector_Free(cones
[i
].vertices
[j
].first
);
450 Matrix_Free(Constraints
);
459 gen_fun
*gen_fun::Hadamard_product(const gen_fun
*gf
, barvinok_options
*options
)
462 Polyhedron
*C
= DomainIntersection(context
, gf
->context
, options
->MaxRays
);
463 Polyhedron
*U
= Universe_Polyhedron(C
->Dimension
);
464 gen_fun
*sum
= new gen_fun(C
);
465 for (short_rat_list::iterator i
= term
.begin(); i
!= term
.end(); ++i
) {
466 for (short_rat_list::iterator i2
= gf
->term
.begin(); i2
!= gf
->term
.end();
468 int d
= (*i
)->d
.power
.NumCols();
469 int k1
= (*i
)->d
.power
.NumRows();
470 int k2
= (*i2
)->d
.power
.NumRows();
471 assert((*i
)->d
.power
.NumCols() == (*i2
)->d
.power
.NumCols());
473 parallel_polytopes
pp((*i
)->n
.power
.NumRows() *
474 (*i2
)->n
.power
.NumRows(),
475 sum
->context
, d
, options
);
477 for (int j
= 0; j
< (*i
)->n
.power
.NumRows(); ++j
) {
478 for (int j2
= 0; j2
< (*i2
)->n
.power
.NumRows(); ++j2
) {
479 Matrix
*M
= Matrix_Alloc(k1
+k2
+d
+d
, 1+k1
+k2
+d
+1);
480 for (int k
= 0; k
< k1
+k2
; ++k
) {
481 value_set_si(M
->p
[k
][0], 1);
482 value_set_si(M
->p
[k
][1+k
], 1);
484 for (int k
= 0; k
< d
; ++k
) {
485 value_set_si(M
->p
[k1
+k2
+k
][1+k1
+k2
+k
], -1);
486 zz2value((*i
)->n
.power
[j
][k
], M
->p
[k1
+k2
+k
][1+k1
+k2
+d
]);
487 for (int l
= 0; l
< k1
; ++l
)
488 zz2value((*i
)->d
.power
[l
][k
], M
->p
[k1
+k2
+k
][1+l
]);
490 for (int k
= 0; k
< d
; ++k
) {
491 value_set_si(M
->p
[k1
+k2
+d
+k
][1+k1
+k2
+k
], -1);
492 zz2value((*i2
)->n
.power
[j2
][k
],
493 M
->p
[k1
+k2
+d
+k
][1+k1
+k2
+d
]);
494 for (int l
= 0; l
< k2
; ++l
)
495 zz2value((*i2
)->d
.power
[l
][k
],
496 M
->p
[k1
+k2
+d
+k
][1+k1
+l
]);
498 Polyhedron
*P
= Constraints2Polyhedron(M
, options
->MaxRays
);
501 QQ c
= (*i
)->n
.coeff
[j
];
502 c
*= (*i2
)->n
.coeff
[j2
];
504 gen_fun
*t
= barvinok_series(P
, U
, options
->MaxRays
);
513 gen_fun
*t
= pp
.compute();
524 void gen_fun::add_union(gen_fun
*gf
, barvinok_options
*options
)
526 QQ
one(1, 1), mone(-1, 1);
528 gen_fun
*hp
= Hadamard_product(gf
, options
);
534 static void Polyhedron_Shift(Polyhedron
*P
, Vector
*offset
)
538 for (int i
= 0; i
< P
->NbConstraints
; ++i
) {
539 Inner_Product(P
->Constraint
[i
]+1, offset
->p
, P
->Dimension
, &tmp
);
540 value_subtract(P
->Constraint
[i
][1+P
->Dimension
],
541 P
->Constraint
[i
][1+P
->Dimension
], tmp
);
543 for (int i
= 0; i
< P
->NbRays
; ++i
) {
544 if (value_notone_p(P
->Ray
[i
][0]))
546 if (value_zero_p(P
->Ray
[i
][1+P
->Dimension
]))
548 Vector_Combine(P
->Ray
[i
]+1, offset
->p
, P
->Ray
[i
]+1,
549 P
->Ray
[i
][0], P
->Ray
[i
][1+P
->Dimension
], P
->Dimension
);
554 void gen_fun::shift(const vec_ZZ
& offset
)
556 for (short_rat_list::iterator i
= term
.begin(); i
!= term
.end(); ++i
)
557 for (int j
= 0; j
< (*i
)->n
.power
.NumRows(); ++j
)
558 (*i
)->n
.power
[j
] += offset
;
560 Vector
*v
= Vector_Alloc(offset
.length());
561 zz2values(offset
, v
->p
);
562 Polyhedron_Shift(context
, v
);
566 /* Divide the generating functin by 1/(1-z^power).
567 * The effect on the corresponding explicit function f(x) is
568 * f'(x) = \sum_{i=0}^\infty f(x - i * power)
570 void gen_fun::divide(const vec_ZZ
& power
)
572 for (short_rat_list::iterator i
= term
.begin(); i
!= term
.end(); ++i
) {
573 int r
= (*i
)->d
.power
.NumRows();
574 int c
= (*i
)->d
.power
.NumCols();
575 (*i
)->d
.power
.SetDims(r
+1, c
);
576 (*i
)->d
.power
[r
] = power
;
579 Vector
*v
= Vector_Alloc(1+power
.length()+1);
580 value_set_si(v
->p
[0], 1);
581 zz2values(power
, v
->p
+1);
582 Polyhedron
*C
= AddRays(v
->p
, 1, context
, context
->NbConstraints
+1);
584 Polyhedron_Free(context
);
588 static void print_power(std::ostream
& os
, QQ
& c
, vec_ZZ
& p
,
589 unsigned int nparam
, char **param_name
)
593 for (int i
= 0; i
< p
.length(); ++i
) {
597 if (c
.n
== -1 && c
.d
== 1)
599 else if (c
.n
!= 1 || c
.d
!= 1) {
615 os
<< "^(" << p
[i
] << ")";
626 void gen_fun::print(std::ostream
& os
, unsigned int nparam
, char **param_name
) const
629 for (short_rat_list::iterator i
= term
.begin(); i
!= term
.end(); ++i
) {
630 if (i
!= term
.begin())
633 for (int j
= 0; j
< (*i
)->n
.coeff
.length(); ++j
) {
634 if (j
!= 0 && (*i
)->n
.coeff
[j
].n
> 0)
636 print_power(os
, (*i
)->n
.coeff
[j
], (*i
)->n
.power
[j
],
640 for (int j
= 0; j
< (*i
)->d
.power
.NumRows(); ++j
) {
644 print_power(os
, mone
, (*i
)->d
.power
[j
], nparam
, param_name
);
651 gen_fun::operator evalue
*() const
655 value_init(factor
.d
);
656 value_init(factor
.x
.n
);
657 for (short_rat_list::iterator i
= term
.begin(); i
!= term
.end(); ++i
) {
658 unsigned nvar
= (*i
)->d
.power
.NumRows();
659 unsigned nparam
= (*i
)->d
.power
.NumCols();
660 Matrix
*C
= Matrix_Alloc(nparam
+ nvar
, 1 + nvar
+ nparam
+ 1);
661 mat_ZZ
& d
= (*i
)->d
.power
;
662 Polyhedron
*U
= context
? context
: Universe_Polyhedron(nparam
);
664 for (int j
= 0; j
< (*i
)->n
.coeff
.length(); ++j
) {
665 for (int r
= 0; r
< nparam
; ++r
) {
666 value_set_si(C
->p
[r
][0], 0);
667 for (int c
= 0; c
< nvar
; ++c
) {
668 zz2value(d
[c
][r
], C
->p
[r
][1+c
]);
670 Vector_Set(&C
->p
[r
][1+nvar
], 0, nparam
);
671 value_set_si(C
->p
[r
][1+nvar
+r
], -1);
672 zz2value((*i
)->n
.power
[j
][r
], C
->p
[r
][1+nvar
+nparam
]);
674 for (int r
= 0; r
< nvar
; ++r
) {
675 value_set_si(C
->p
[nparam
+r
][0], 1);
676 Vector_Set(&C
->p
[nparam
+r
][1], 0, nvar
+ nparam
+ 1);
677 value_set_si(C
->p
[nparam
+r
][1+r
], 1);
679 Polyhedron
*P
= Constraints2Polyhedron(C
, 0);
680 evalue
*E
= barvinok_enumerate_ev(P
, U
, 0);
682 if (EVALUE_IS_ZERO(*E
)) {
687 zz2value((*i
)->n
.coeff
[j
].n
, factor
.x
.n
);
688 zz2value((*i
)->n
.coeff
[j
].d
, factor
.d
);
691 Matrix_Print(stdout, P_VALUE_FMT, C);
692 char *test[] = { "A", "B", "C", "D", "E", "F", "G" };
693 print_evalue(stdout, E, test);
707 value_clear(factor
.d
);
708 value_clear(factor
.x
.n
);
712 void gen_fun::coefficient(Value
* params
, Value
* c
) const
714 if (context
&& !in_domain(context
, params
)) {
721 value_init(part
.x
.n
);
724 evalue_set_si(&sum
, 0, 1);
728 for (short_rat_list::iterator i
= term
.begin(); i
!= term
.end(); ++i
) {
729 unsigned nvar
= (*i
)->d
.power
.NumRows();
730 unsigned nparam
= (*i
)->d
.power
.NumCols();
731 Matrix
*C
= Matrix_Alloc(nparam
+ nvar
, 1 + nvar
+ 1);
732 mat_ZZ
& d
= (*i
)->d
.power
;
734 for (int j
= 0; j
< (*i
)->n
.coeff
.length(); ++j
) {
735 C
->NbRows
= nparam
+nvar
;
736 for (int r
= 0; r
< nparam
; ++r
) {
737 value_set_si(C
->p
[r
][0], 0);
738 for (int c
= 0; c
< nvar
; ++c
) {
739 zz2value(d
[c
][r
], C
->p
[r
][1+c
]);
741 zz2value((*i
)->n
.power
[j
][r
], C
->p
[r
][1+nvar
]);
742 value_subtract(C
->p
[r
][1+nvar
], C
->p
[r
][1+nvar
], params
[r
]);
744 for (int r
= 0; r
< nvar
; ++r
) {
745 value_set_si(C
->p
[nparam
+r
][0], 1);
746 Vector_Set(&C
->p
[nparam
+r
][1], 0, nvar
+ 1);
747 value_set_si(C
->p
[nparam
+r
][1+r
], 1);
749 Polyhedron
*P
= Constraints2Polyhedron(C
, 0);
754 barvinok_count(P
, &tmp
, 0);
756 if (value_zero_p(tmp
))
758 zz2value((*i
)->n
.coeff
[j
].n
, part
.x
.n
);
759 zz2value((*i
)->n
.coeff
[j
].d
, part
.d
);
760 value_multiply(part
.x
.n
, part
.x
.n
, tmp
);
766 assert(value_one_p(sum
.d
));
767 value_assign(*c
, sum
.x
.n
);
771 value_clear(part
.x
.n
);
773 value_clear(sum
.x
.n
);
776 gen_fun
*gen_fun::summate(int nvar
, barvinok_options
*options
) const
778 int dim
= context
->Dimension
;
779 int nparam
= dim
- nvar
;
783 if (options
->incremental_specialization
== 1) {
784 red
= new partial_ireducer(Polyhedron_Project(context
, nparam
), dim
, nparam
);
786 red
= new partial_reducer(Polyhedron_Project(context
, nparam
), dim
, nparam
);
788 for (short_rat_list::iterator i
= term
.begin(); i
!= term
.end(); ++i
)
789 for (int j
= 0; j
< (*i
)->n
.power
.NumRows(); ++j
)
790 red
->reduce((*i
)->n
.coeff
[j
], (*i
)->n
.power
[j
], (*i
)->d
.power
);
796 /* returns true if the set was finite and false otherwise */
797 bool gen_fun::summate(Value
*sum
) const
799 if (term
.size() == 0) {
800 value_set_si(*sum
, 0);
805 for (short_rat_list::iterator i
= term
.begin(); i
!= term
.end(); ++i
)
806 if ((*i
)->d
.power
.NumRows() > maxlen
)
807 maxlen
= (*i
)->d
.power
.NumRows();
809 infinite_icounter
cnt((*term
.begin())->d
.power
.NumCols(), maxlen
);
810 for (short_rat_list::iterator i
= term
.begin(); i
!= term
.end(); ++i
)
811 for (int j
= 0; j
< (*i
)->n
.power
.NumRows(); ++j
)
812 cnt
.reduce((*i
)->n
.coeff
[j
], (*i
)->n
.power
[j
], (*i
)->d
.power
);
814 for (int i
= 1; i
<= maxlen
; ++i
)
815 if (value_notzero_p(mpq_numref(cnt
.count
[i
]))) {
816 value_set_si(*sum
, -1);
820 assert(value_one_p(mpq_denref(cnt
.count
[0])));
821 value_assign(*sum
, mpq_numref(cnt
.count
[0]));