2 #include <piplib/piplibMP.h>
3 #include <barvinok/basis_reduction.h>
13 static struct pip_lp
*init_lp(Polyhedron
*P
);
14 static void set_lp_obj(struct pip_lp
*lp
, Value
*row
, int dim
);
15 static int solve_lp(struct pip_lp
*lp
);
16 static void get_obj_val(struct pip_lp
* lp
, mpq_t
*F
);
17 static void delete_lp(struct pip_lp
*lp
);
18 static int add_lp_row(struct pip_lp
*lp
, Value
*row
, int dim
);
19 static void get_alpha(struct pip_lp
* lp
, int row
, mpq_t
*alpha
);
20 static void del_lp_row(struct pip_lp
*lp
);
22 #define GBR_LP struct pip_lp
23 #define GBR_type mpq_t
24 #define GBR_init(v) mpq_init(v)
25 #define GBR_clear(v) mpq_clear(v)
26 #define GBR_set(a,b) mpq_set(a,b)
27 #define GBR_set_ui(a,b) mpq_set_ui(a,b,1)
28 #define GBR_mul(a,b,c) mpq_mul(a,b,c)
29 #define GBR_lt(a,b) (mpq_cmp(a,b) < 0)
30 #define GBR_floor(a,b) mpz_fdiv_q(a,mpq_numref(b),mpq_denref(b))
31 #define GBR_ceil(a,b) mpz_cdiv_q(a,mpq_numref(b),mpq_denref(b))
32 #define GBR_lp_init(P) init_lp(P)
33 #define GBR_lp_set_obj(lp, obj, dim) set_lp_obj(lp, obj, dim)
34 #define GBR_lp_solve(lp) solve_lp(lp)
35 #define GBR_lp_get_obj_val(lp, F) get_obj_val(lp, F)
36 #define GBR_lp_delete(lp) delete_lp(lp)
37 #define GBR_lp_next_row(lp) lp->neq
38 #define GBR_lp_add_row(lp, row, dim) add_lp_row(lp, row, dim)
39 #define GBR_lp_get_alpha(lp, row, alpha) get_alpha(lp, row, alpha)
40 #define GBR_lp_del_row(lp) del_lp_row(lp);
41 #define Polyhedron_Reduced_Basis pip_dual_Polyhedron_Reduced_Basis
42 #include "basis_reduction_templ.c"
44 #define ALLOC(type) (type*)malloc(sizeof(type))
46 static struct pip_lp
*init_lp(Polyhedron
*P
)
48 struct pip_lp
*lp
= ALLOC(struct pip_lp
);
50 lp
->eq
= Matrix_Alloc(P
->Dimension
, P
->Dimension
);
57 static void set_lp_obj(struct pip_lp
*lp
, Value
*row
, int dim
)
59 assert(lp
->P
->Dimension
== dim
);
63 static int solve_lp(struct pip_lp
*lp
)
67 unsigned dim
= lp
->P
->Dimension
;
68 int rows
= 1 + 2*dim
+ 2*lp
->P
->NbConstraints
;
69 int cols
= 1 + 1 + lp
->neq
+ 2*lp
->P
->NbConstraints
+ 1;
73 pip_quast_free(lp
->sol
);
75 assert(lp
->P
->NbEq
== 0);
76 domain
= pip_matrix_alloc(rows
, cols
);
77 for (i
= 0; i
< lp
->neq
; ++i
)
78 for (j
= 0; j
< dim
; ++j
) {
79 value_assign(domain
->p
[j
][1+1+i
], lp
->eq
->p
[i
][j
]);
80 value_oppose(domain
->p
[dim
+j
][1+1+i
], lp
->eq
->p
[i
][j
]);
82 value_set_si(domain
->p
[2*dim
][1], 1);
83 for (j
= 0; j
< dim
; ++j
) {
84 value_oppose(domain
->p
[j
][cols
-1], lp
->obj
[j
]);
85 value_assign(domain
->p
[dim
+j
][cols
-1], lp
->obj
[j
]);
87 for (i
= 0; i
< lp
->P
->NbConstraints
; ++i
) {
88 for (j
= 0; j
< dim
; ++j
) {
89 value_assign(domain
->p
[j
][1+1+lp
->neq
+i
], lp
->P
->Constraint
[i
][1+j
]);
90 value_assign(domain
->p
[dim
+j
][1+1+lp
->neq
+lp
->P
->NbConstraints
+i
],
91 lp
->P
->Constraint
[i
][1+j
]);
93 value_assign(domain
->p
[2*dim
][1+1+lp
->neq
+i
], lp
->P
->Constraint
[i
][1+dim
]);
94 value_assign(domain
->p
[2*dim
][1+1+lp
->neq
+lp
->P
->NbConstraints
+i
],
95 lp
->P
->Constraint
[i
][1+dim
]);
96 value_set_si(domain
->p
[2*dim
+1+i
][0], 1);
97 value_set_si(domain
->p
[2*dim
+1+i
][1+1+lp
->neq
+i
], -1);
98 value_set_si(domain
->p
[2*dim
+1+lp
->P
->NbConstraints
+i
][0], 1);
99 value_set_si(domain
->p
[2*dim
+1+lp
->P
->NbConstraints
+i
]
100 [1+1+lp
->neq
+lp
->P
->NbConstraints
+i
], -1);
103 options
= pip_options_init();
104 options
->Urs_unknowns
= -1;
106 lp
->sol
= pip_solve(domain
, NULL
, -1, options
);
107 pip_matrix_free(domain
);
109 pip_options_free(options
);
111 return value_zero_p(lp
->sol
->list
->vector
->the_deno
[0]);
114 static void get_obj_val(struct pip_lp
* lp
, mpq_t
*F
)
116 value_assign(mpq_numref(*F
), lp
->sol
->list
->vector
->the_vector
[0]);
117 value_assign(mpq_denref(*F
), lp
->sol
->list
->vector
->the_deno
[0]);
120 static void delete_lp(struct pip_lp
*lp
)
123 pip_quast_free(lp
->sol
);
128 static int add_lp_row(struct pip_lp
*lp
, Value
*row
, int dim
)
130 assert(lp
->P
->Dimension
== dim
);
131 Vector_Copy(row
, lp
->eq
->p
[lp
->neq
], dim
);
135 static void get_alpha(struct pip_lp
* lp
, int row
, mpq_t
*alpha
)
140 for (i
= 0, list
= lp
->sol
->list
->next
; i
< row
; ++i
, list
= list
->next
)
143 if (value_zero_p(list
->vector
->the_deno
[0])) {
144 value_set_si(mpq_numref(*alpha
), 0);
145 value_set_si(mpq_denref(*alpha
), 1);
147 value_oppose(mpq_numref(*alpha
), list
->vector
->the_vector
[0]);
148 value_assign(mpq_denref(*alpha
), list
->vector
->the_deno
[0]);
152 static void del_lp_row(struct pip_lp
*lp
)