1 #include <polylib/polylibgmp.h>
2 #include <stdio.h> /* needed for piplib ! */
3 #include <piplib/piplibMP.h>
8 #ifdef HAVE_GROWING_CHERNIKOVA
9 #define MAXRAYS POL_NO_DUAL
14 static PipMatrix
*poly2pip(Polyhedron
*P
, int pos
, int n
, int nparam
)
18 unsigned int extra
= P
->Dimension
- pos
- n
- nparam
;
20 matrix
= pip_matrix_alloc(P
->NbConstraints
, P
->Dimension
+2);
21 for (i
= 0; i
< P
->NbConstraints
; ++i
) {
22 value_assign(matrix
->p
[i
][0], P
->Constraint
[i
][0]);
23 for (j
= 0; j
< n
; ++j
)
24 value_assign(matrix
->p
[i
][1+j
], P
->Constraint
[i
][1+pos
+j
]);
25 for (j
= 0; j
< pos
; ++j
)
26 value_assign(matrix
->p
[i
][1+n
+j
], P
->Constraint
[i
][1+j
]);
27 for (j
= 0; j
< extra
; ++j
)
28 value_assign(matrix
->p
[i
][1+n
+pos
+j
], P
->Constraint
[i
][1+n
+pos
+j
]);
29 for (j
= 1+pos
+n
+extra
; j
< P
->Dimension
+2; ++j
)
30 value_assign(matrix
->p
[i
][j
], P
->Constraint
[i
][j
]);
36 static int max_new(PipQuast
*q
, int max
, int d
, int *maxd
)
40 for (p
= q
->newparm
; p
; p
= p
->next
)
46 max
= max_new(q
->next_else
, max
, d
, maxd
);
47 max
= max_new(q
->next_then
, max
, d
, maxd
);
53 * v: total number of variables in original problem
54 * pos: position of minimized variables
55 * n: dimension of space in which minimum was taken
56 * i.e., number of minimized vars in original problem
57 * e: number of extra existential variables
67 void (*f
)(struct scan_data
*data
, PipQuast
*q
, int i
);
70 static int vectorpos2cpos(struct scan_data
*data
, int j
)
72 int np
= data
->d
- data
->v
- data
->e
;
77 else if (j
< data
->v
- data
->n
)
79 else if (j
< data
->v
-data
->n
+np
)
80 l
= 1+data
->n
+data
->e
+j
;
90 * Dimensions in the resulting domain are:
93 static void scan_quast_r(struct scan_data
*data
, PipQuast
*q
, int i
)
97 for (p
= q
->newparm
; p
; p
= p
->next
) {
99 PipVector
*vec
= p
->vector
;
101 value_set_si(data
->M
->p
[i
][0], 1);
102 Vector_Set(data
->M
->p
[i
]+1, 0, data
->d
);
104 for (j
= 0; j
< vec
->nb_elements
-1; ++j
) {
105 l
= vectorpos2cpos(data
, j
);
106 value_assign(data
->M
->p
[i
][l
], vec
->the_vector
[j
]);
108 l
= vectorpos2cpos(data
, p
->rank
);
109 value_oppose(data
->M
->p
[i
][l
], p
->deno
);
110 value_assign(data
->M
->p
[i
][1+data
->d
], vec
->the_vector
[j
]);
113 value_set_si(data
->M
->p
[i
][0], 1);
114 Vector_Set(data
->M
->p
[i
]+1, 0, data
->d
);
116 for (j
= 0; j
< vec
->nb_elements
-1; ++j
) {
117 l
= vectorpos2cpos(data
, j
);
118 value_oppose(data
->M
->p
[i
][l
], vec
->the_vector
[j
]);
120 l
= vectorpos2cpos(data
, p
->rank
);
121 value_assign(data
->M
->p
[i
][l
], p
->deno
);
122 value_oppose(data
->M
->p
[i
][1+data
->d
], vec
->the_vector
[j
]);
123 value_addto(data
->M
->p
[i
][1+data
->d
], data
->M
->p
[i
][1+data
->d
], p
->deno
);
124 value_decrement(data
->M
->p
[i
][1+data
->d
], data
->M
->p
[i
][1+data
->d
]);
130 value_set_si(data
->M
->p
[i
][0], 1);
131 Vector_Set(data
->M
->p
[i
]+1, 0, data
->d
);
133 for (j
= 0; j
< q
->condition
->nb_elements
-1; ++j
) {
134 int l
= vectorpos2cpos(data
, j
);
135 value_assign(data
->M
->p
[i
][l
], q
->condition
->the_vector
[j
]);
137 value_assign(data
->M
->p
[i
][1+data
->d
], q
->condition
->the_vector
[j
]);
138 scan_quast_r(data
, q
->next_then
, i
+1);
140 for (j
= 0; j
< q
->condition
->nb_elements
-1; ++j
) {
141 int l
= vectorpos2cpos(data
, j
);
142 value_oppose(data
->M
->p
[i
][l
], q
->condition
->the_vector
[j
]);
144 value_oppose(data
->M
->p
[i
][1+data
->d
], q
->condition
->the_vector
[j
]);
145 value_decrement(data
->M
->p
[i
][1+data
->d
], data
->M
->p
[i
][1+data
->d
]);
146 scan_quast_r(data
, q
->next_else
, i
+1);
151 /* if data->n is zero, we are only interested in the domains
152 * where a solution exists and not in the actual solution
154 if (q
->list
&& data
->n
) {
157 for (j
= 0, l
= q
->list
; l
; ++j
, l
= l
->next
) {
158 Vector_Set(data
->M
->p
[i
], 0, data
->d
+1);
159 value_set_si(data
->M
->p
[i
][1+data
->pos
+j
], -1);
161 for (k
= 0; k
< l
->vector
->nb_elements
-1; ++k
) {
162 int ll
= vectorpos2cpos(data
, k
);
163 value_assign(data
->M
->p
[i
][ll
], l
->vector
->the_vector
[k
]);
165 value_assign(data
->M
->p
[i
][1+data
->d
], l
->vector
->the_vector
[k
]);
173 static void scan_quast(struct scan_data
*data
, PipQuast
*q
)
175 int nparam
, nexist
, d
, i
, rows
;
177 nparam
= data
->d
- data
->n
;
180 nexist
= max_new(q
, nparam
-1, 0, &d
) - nparam
+1;
181 rows
= nparam
+ 2 * nexist
+ d
+ data
->n
;
183 /* nparam now refers to the number of parameters in the original polyhedron */
184 nparam
-= data
->v
- data
->n
;
187 data
->d
= data
->v
+ nexist
+ nparam
;
189 data
->M
= Matrix_Alloc(rows
, 1+data
->d
+1);
190 /* All parameters are/should be positive */
191 for (i
= 0; i
< data
->pos
; ++i
) {
192 value_set_si(data
->M
->p
[i
][0], 1);
193 value_set_si(data
->M
->p
[i
][1+i
], 1);
195 for (i
= data
->pos
+data
->n
; i
< data
->v
; ++i
) {
196 value_set_si(data
->M
->p
[i
-data
->n
][0], 1);
197 value_set_si(data
->M
->p
[i
-data
->n
][1+i
], 1);
199 for (i
= 0; i
< nparam
; ++i
) {
200 value_set_si(data
->M
->p
[data
->v
-data
->n
+i
][0], 1);
201 value_set_si(data
->M
->p
[data
->v
-data
->n
+i
][1+data
->v
+nexist
+i
], 1);
204 scan_quast_r(data
, q
, data
->d
- data
->e
- data
->n
);
205 Matrix_Free(data
->M
);
208 struct quast2poly_data
{
209 struct scan_data scan
;
213 static void add_quast_base(struct scan_data
*data
, PipQuast
*q
, int i
)
215 struct quast2poly_data
*qdata
= (struct quast2poly_data
*)data
;
219 C
= Matrix_Alloc(i
, 1+data
->d
+1);
220 Vector_Copy(data
->M
->p
[0], C
->p
[0], i
* (1+data
->d
+1));
221 P
= Constraints2Polyhedron(C
, MAXRAYS
);
229 * nvar: total number of variables, including the minimized variables,
230 * but excluding the parameters
231 * nparam: the number of parameters in the PIP problem, excluding the "big parameter"
232 * for maximization problems, i.e., the
233 * total number of variables minus the minimized variables
234 * pos: position of the first minimized variable
235 * n: number of minimized variables
237 static Polyhedron
*quast2poly(PipQuast
*q
, int nvar
, int nparam
, int pos
, int n
)
239 struct quast2poly_data data
;
244 data
.scan
.d
= n
+nparam
;
247 data
.scan
.f
= add_quast_base
;
249 scan_quast(&data
.scan
, q
);
254 Polyhedron
*pip_projectout(Polyhedron
*P
, int pos
, int n
, int nparam
)
258 unsigned int nvar
= P
->Dimension
- nparam
;
259 PipMatrix
*context
= pip_matrix_alloc(0, P
->Dimension
- n
+ 2);
263 POL_ENSURE_INEQUALITIES(P
);
265 domain
= poly2pip(P
, pos
, n
, nparam
);
267 pip_matrix_print(stderr
, domain
);
268 pip_matrix_print(stderr
, context
);
271 options
= pip_options_init();
272 options
->Simplify
= 1;
273 sol
= pip_solve(domain
, context
, -1, options
);
276 pip_quast_print(stderr
, sol
, 0);
279 min
= quast2poly(sol
, nvar
- n
, P
->Dimension
- n
, 0, 0);
282 pip_matrix_free(context
);
283 pip_matrix_free(domain
);
284 pip_options_free(options
);