occ: add interface to barvinok_summate
[barvinok.git] / zsolve / valuetrees.c
blobbdf5d4ed9b6c70fe5829a6617df23e9a4422eea6
1 /*
2 4ti2 -- A software package for algebraic, geometric and combinatorial
3 problems on linear spaces.
5 Copyright (C) 2006 4ti2 team.
6 Main author(s): Matthias Walter.
8 This program is free software; you can redistribute it and/or
9 modify it under the terms of the GNU General Public License
10 as published by the Free Software Foundation; either version 2
11 of the License, or (at your option) any later version.
13 This program is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
18 You should have received a copy of the GNU General Public License
19 along with this program; if not, write to the Free Software
20 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
23 #include "defs.h"
24 #include "vectorarray.h"
25 #include "indexarray.h"
26 #include "valuetrees.h"
28 #include <assert.h>
29 #include <stdlib.h>
30 #include <stdio.h>
32 // //
34 ValueTree createValueTree(int level)
36 ValueTree tree = (ValueTree)malloc(sizeof(valuetree_t));
38 if (tree==NULL)
40 fprintf(stderr, "Fatal Error: Could not allocate memory for ValueTree!\n");
41 exit(1);
44 tree->level = level;
45 tree->pos = NULL;
46 tree->neg = NULL;
47 tree->zero = NULL;
48 tree->vectors = level<0 ? createIndexArray() : NULL;
50 return tree;
53 // //
55 void deleteValueTree(ValueTree tree)
57 ValueTreeNode node;
59 if (tree!=NULL)
61 while ((node = tree->pos)!=NULL)
63 tree->pos = node->next;
64 deleteValueTree(node->sub);
65 free(node);
67 while ((node = tree->neg)!=NULL)
69 tree->neg = node->next;
70 deleteValueTree(node->sub);
71 free(node);
73 deleteValueTree(tree->zero);
74 deleteIndexArray(tree->vectors);
75 free(tree);
79 // //
81 void splitValueTree(ZSolveContext ctx, ValueTree tree, int start)
83 int i, compo, value;
84 bool pos, neg, zero;
85 ValueTreeNode node, temp;
87 if (!tree || tree->level>=0)
88 return;
90 assert(tree->pos==NULL);
91 assert(tree->neg==NULL);
92 assert(tree->zero==NULL);
93 assert(tree->vectors);
96 for (; start<ctx->Current; start++)
98 compo = start<0 ? ctx->Current : start;
99 pos = false;
100 neg = false;
101 zero = false;
102 for (i=0; i<tree->vectors->Size; i++)
104 if (ctx->Lattice->Data[tree->vectors->Data[i]][compo]>0)
105 pos = true;
106 else if (ctx->Lattice->Data[tree->vectors->Data[i]][compo]<0)
107 neg = true;
108 else
109 zero = true;
111 if (pos+neg>1)
112 break;
115 if (start<ctx->Current && tree->vectors->Size > 1)
117 compo = start<0 ? ctx->Current : start;
118 for (i=0; i<tree->vectors->Size; i++)
120 value = ctx->Lattice->Data[tree->vectors->Data[i]][compo];
121 if (value>0)
123 node = tree->pos;
124 if (node==NULL || value<node->value)
126 tree->pos = (ValueTreeNode)malloc(sizeof(valuetreenode_t));
127 tree->pos->value = value;
128 tree->pos->next = node;
129 tree->pos->sub = createValueTree(-1);
130 appendToIndexArray(tree->pos->sub->vectors, tree->vectors->Data[i]);
132 else
134 while (node->next && value>=node->next->value)
135 node = node->next;
136 if (value!=node->value)
138 temp = (ValueTreeNode)malloc(sizeof(valuetreenode_t));
139 temp->value = value;
140 temp->sub = createValueTree(-1);
141 temp->next = node->next;
142 node->next = temp;
143 node = temp;
145 appendToIndexArray(node->sub->vectors, tree->vectors->Data[i]);
148 else if (value<0)
150 node = tree->neg;
151 if (node==NULL || value>node->value)
153 tree->neg = (ValueTreeNode)malloc(sizeof(valuetreenode_t));
154 tree->neg->value = value;
155 tree->neg->next = node;
156 tree->neg->sub = createValueTree(-1);
157 appendToIndexArray(tree->neg->sub->vectors, tree->vectors->Data[i]);
159 else
161 while (node->next && value<=node->next->value)
162 node = node->next;
163 if (value!=node->value)
165 temp = (ValueTreeNode)malloc(sizeof(valuetreenode_t));
166 temp->value = value;
167 temp->sub = createValueTree(-1);
168 temp->next = node->next;
169 node->next = temp;
170 node = temp;
172 appendToIndexArray(node->sub->vectors, tree->vectors->Data[i]);
175 else
177 if (tree->zero==NULL)
178 tree->zero = createValueTree(-1);
179 appendToIndexArray(tree->zero->vectors, tree->vectors->Data[i]);
183 deleteIndexArray(tree->vectors);
184 tree->vectors = NULL;
185 tree->level = compo;
186 start++;
187 node = tree->pos;
188 while (node)
190 splitValueTree(ctx, node->sub, start);
191 node = node->next;
193 node = tree->neg;
194 while (node)
196 splitValueTree(ctx, node->sub, start);
197 node = node->next;
200 splitValueTree(ctx, tree->zero, start);
201 return;
205 // //
207 void createValueTrees(ZSolveContext ctx)
209 int i,j;
211 ctx->MaxNorm = 0;
212 for (i=0; i<ctx->Lattice->Size; i++)
213 ctx->MaxNorm = maxi(ctx->MaxNorm, normVector(ctx->Lattice->Data[i], ctx->Current));
215 ctx->Norm = (void **)calloc(ctx->MaxNorm+1, sizeof(ValueTree));
217 for (i=0; i<ctx->Lattice->Size; i++)
219 j = normVector(ctx->Lattice->Data[i], ctx->Current);
220 // TODO is it correct, to throw out if all 0 on before and current ?! - it seems incorrect, as bayer-test fails.
221 if (j==0 && ctx->Lattice->Data[i][ctx->Current] == 0)
222 continue;
223 // put in corresponding norm-tree
224 if (ctx->Norm[j]==NULL)
225 ctx->Norm[j] = createValueTree(-1);
226 appendToIndexArray( ((ValueTree *)ctx->Norm)[j]->vectors, i);
229 for (i=0; i<=ctx->MaxNorm; i++)
230 splitValueTree(ctx, ctx->Norm[i], -1);
233 // //
235 void deleteValueTrees(ZSolveContext ctx)
237 int i;
239 for (i=0; i<=ctx->MaxNorm; i++)
240 deleteValueTree(ctx->Norm[i]);
242 free(ctx->Norm);
243 ctx->Norm = NULL;
246 // //
248 bool enumValueReducer(ZSolveContext ctx, ValueTree tree)
250 int value, i, j;
251 Vector Reducer;
252 ValueTreeNode node;
254 if (tree)
256 if (tree->level>=0)
258 value = ctx->Sum[tree->level];
259 if (value>0)
261 node = tree->pos;
262 while (node && node->value<=value)
264 if (enumValueReducer(ctx, node->sub))
265 return true;
266 node = node->next;
269 if (value<0)
271 node = tree->neg;
272 while (node && node->value>=value)
274 if (enumValueReducer(ctx, node->sub))
275 return true;
276 node = node->next;
279 if (enumValueReducer(ctx, tree->zero))
280 return true;
282 else
284 for (i=0; i<tree->vectors->Size; i++)
286 Reducer = ctx->Lattice->Data[tree->vectors->Data[i]];
287 /* printf("=> Testing reduction of [");
288 printVector(ctx->Sum, ctx->Variables);
289 printf("] by [");
290 printVector(Reducer, ctx->Variables);
291 printf("]\n"); */
292 for (j=0; j<=ctx->Current; j++)
293 if (Reducer[j]*ctx->Sum[j]<0 || abs(ctx->Sum[j])<abs(Reducer[j]))
294 break;
295 if (j>ctx->Current)
296 return true;
300 return false;
303 // //
305 void insertVectorToValueTree(ZSolveContext ctx, ValueTree tree, int index)
307 int value;
308 ValueTreeNode node, temp;
310 assert(tree);
311 assert(index>=0 && index<ctx->Lattice->Size);
313 if (tree->level>=0)
315 value = ctx->Lattice->Data[index][tree->level];
316 if (value>0)
318 node = tree->pos;
319 if (node==NULL || value<node->value)
321 tree->pos = (ValueTreeNode)malloc(sizeof(valuetreenode_t));
322 tree->pos->value = value;
323 tree->pos->next = node;
324 tree->pos->sub = createValueTree(-1);
325 appendToIndexArray(tree->pos->sub->vectors, index);
327 else
329 while (node->next && value>=node->next->value)
330 node = node->next;
331 if (value==node->value)
332 insertVectorToValueTree(ctx, node->sub, index);
333 else
335 temp = (ValueTreeNode)malloc(sizeof(valuetreenode_t));
336 temp->value = value;
337 temp->sub = createValueTree(-1);
338 appendToIndexArray(temp->sub->vectors, index);
339 temp->next = node->next;
340 node->next = temp;
344 else if (value<0)
346 node = tree->neg;
347 if (node==NULL || value>node->value)
349 tree->neg = (ValueTreeNode)malloc(sizeof(valuetreenode_t));
350 tree->neg->value = value;
351 tree->neg->next = node;
352 tree->neg->sub = createValueTree(-1);
353 appendToIndexArray(tree->neg->sub->vectors, index);
355 else
357 while (node->next && value<=node->next->value)
358 node = node->next;
359 if (value==node->value)
360 insertVectorToValueTree(ctx, node->sub, index);
361 else
363 temp = (ValueTreeNode)malloc(sizeof(valuetreenode_t));
364 temp->value = value;
365 temp->sub = createValueTree(-1);
366 appendToIndexArray(temp->sub->vectors, index);
367 temp->next = node->next;
368 node->next = temp;
372 else
374 if (tree->zero==NULL)
376 tree->zero = createValueTree(-1);
377 appendToIndexArray(tree->zero->vectors, index);
379 else
380 insertVectorToValueTree(ctx, tree->zero, index);
383 else
385 appendToIndexArray(tree->vectors, index);
386 splitValueTree(ctx, tree, -1);
390 // //
392 void insertVectorToValueTrees(ZSolveContext ctx, Vector vector, int norm)
394 if (norm>ctx->MaxNorm)
396 ctx->Norm = (void **)realloc(ctx->Norm, (norm+1)*sizeof(ValueTree));
397 while (norm>ctx->MaxNorm)
398 ctx->Norm[++ctx->MaxNorm] = NULL;
401 appendToVectorArray(ctx->Lattice, vector);
402 if (ctx->Norm[norm]==NULL)
404 ctx->Norm[norm] = createValueTree(-1);
405 appendToIndexArray(((ValueTree *)ctx->Norm)[norm]->vectors, ctx->Lattice->Size-1);
407 else
408 insertVectorToValueTree(ctx, ctx->Norm[norm], ctx->Lattice->Size-1);
411 // //
413 void buildValueSum(ZSolveContext ctx)
415 int i,norm;
416 bool flag;
418 assert(ctx->First);
419 assert(ctx->Second);
420 assert(ctx->Sum);
422 // same ?
423 if (ctx->First==ctx->Second)
424 return;
426 // pattern okay ?
427 if (ctx->First[ctx->Current]*ctx->Second[ctx->Current]>0)
428 return;
429 flag = false;
430 for (i=0; i<ctx->Current; i++)
431 if (ctx->First[i]*ctx->Second[i]<0)
432 flag = true;
433 if (flag)
434 return;
436 // build new
437 norm = 0;
438 for (i=0; i<ctx->Variables; i++)
440 ctx->Sum[i] = ctx->First[i] + ctx->Second[i];
441 if (i<ctx->Current)
442 norm += abs(ctx->Sum[i]);
444 // all zero?
445 if (norm==0)
446 return;
448 // DEBUG
449 /* printf("build: [");
450 printVector(ctx->Sum, ctx->Variables);
451 printf("]\n"); */
453 // reducable ?
454 flag = false;
455 for (i=0; i<=norm/2; i++)
457 if (i<=ctx->MaxNorm)
459 if (enumValueReducer(ctx, ctx->Norm[i]))
461 flag = true;
462 break;
466 if (flag) {
467 // printf("(reduced)\n");
468 return;
471 // limits okay ?
472 for (i=0; i<ctx->Current; i++)
474 if (ctx->Lattice->Properties[i].Lower>ctx->Sum[i]) {
475 // printf("(lower bounds failed)\n");
476 return;
478 if (ctx->Lattice->Properties[i].Upper<ctx->Sum[i]) {
479 // printf("(upper bounds failed)\n");
480 return;
484 if (norm<=ctx->MaxNorm)
485 if (enumValueReducer(ctx, ctx->Norm[norm]))
486 return;
488 insertVectorToValueTrees(ctx, copyVector(ctx->Sum, ctx->Variables), norm);
490 if (ctx->Symmetric)
492 for (i=0; i<ctx->Variables; i++)
493 ctx->Sum[i] = -ctx->Sum[i];
494 insertVectorToValueTrees(ctx, copyVector(ctx->Sum, ctx->Variables), norm);
498 // //
500 void enumValueSecond(ZSolveContext ctx, ValueTree tree)
502 int i;
503 ValueTreeNode node;
505 if (tree)
507 if (tree->level==ctx->Current)
509 if (ctx->First[ctx->Current]>=0)
511 node = tree->neg;
512 while (node)
514 enumValueSecond(ctx, node->sub);
515 node = node->next;
518 if (ctx->First[ctx->Current]<=0)
520 node = tree->pos;
521 while (node)
523 enumValueSecond(ctx, node->sub);
524 node = node->next;
528 else if (tree->level>=0)
530 if (ctx->First[tree->level]<=0)
532 node = tree->neg;
533 while (node)
535 enumValueSecond(ctx, node->sub);
536 node = node->next;
539 enumValueSecond(ctx, tree->zero);
540 if (ctx->First[tree->level]>=0)
542 node = tree->pos;
543 while (node)
545 enumValueSecond(ctx, node->sub);
546 node = node->next;
550 else
552 assert(tree->vectors);
554 for (i=0; tree->vectors!=NULL && i<tree->vectors->Size; i++)
556 ctx->Second = ctx->Lattice->Data[tree->vectors->Data[i]];
557 if (ctx->Second[ctx->Current]!=0)
558 buildValueSum(ctx);
560 if (tree->vectors==NULL)
561 enumValueSecond(ctx, tree);
566 // //
568 void enumValueFirst(ZSolveContext ctx, ValueTree tree, int norm)
570 int i;
572 ValueTreeNode node;
575 if (tree)
577 if (tree->level>=0)
579 node = tree->pos;
580 while (node)
582 enumValueFirst(ctx, node->sub, norm);
583 node = node->next;
585 enumValueFirst(ctx, tree->zero, norm);
586 node = tree->neg;
587 while (node)
589 enumValueFirst(ctx, node->sub, norm);
590 node = node->next;
593 else
595 assert(tree->vectors);
597 for (i=0; tree->vectors!=NULL && i<tree->vectors->Size; i++)
599 ctx->First = ctx->Lattice->Data[tree->vectors->Data[i]];
600 if ((!ctx->Symmetric && ctx->First[ctx->Current]<0) || ctx->First[ctx->Current]>0)
601 enumValueSecond(ctx, ctx->Norm[norm]);
604 if (tree->vectors==NULL)
605 enumValueFirst(ctx, tree, norm);
610 // //
612 void completeValueTrees(ZSolveContext ctx, int norm1, int norm2)
614 assert(ctx->Norm);
615 assert(!(norm1>ctx->MaxNorm || norm2>ctx->MaxNorm || ctx->Norm[norm1]==NULL || ctx->Norm[norm2]==NULL));
617 enumValueFirst(ctx, ctx->Norm[norm1], norm2);
620 // //