From c8ece89973f081b38ee39cdc68481153af6a11a2 Mon Sep 17 00:00:00 2001 From: Sven Verdoolaege Date: Tue, 1 Apr 2008 22:12:56 +0200 Subject: [PATCH] normalization.c: extract standard_constraints from polysign.c --- Makefile.am | 2 ++ hilbert.c | 1 + normalization.c | 88 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++ normalization.h | 12 ++++++++ polysign.c | 86 ------------------------------------------------------- polysign.h | 3 -- topcom.c | 2 +- 7 files changed, 104 insertions(+), 90 deletions(-) create mode 100644 normalization.c create mode 100644 normalization.h diff --git a/Makefile.am b/Makefile.am index 604de58..d997fc7 100644 --- a/Makefile.am +++ b/Makefile.am @@ -113,6 +113,8 @@ libbarvinok_core_la_SOURCES = \ lattice_point.h \ lattice_width.c \ lattice_width.h \ + normalization.c \ + normalization.h \ options.c \ param_util.c \ param_util.h \ diff --git a/hilbert.c b/hilbert.c index 3ea20bd..546e258 100644 --- a/hilbert.c +++ b/hilbert.c @@ -8,6 +8,7 @@ #include #include #include "hilbert.h" +#include "normalization.h" #include "polysign.h" static ZSolveMatrix Matrix2zsolve(Matrix *M) diff --git a/normalization.c b/normalization.c new file mode 100644 index 0000000..d5e9080 --- /dev/null +++ b/normalization.c @@ -0,0 +1,88 @@ +#include +#include "normalization.h" + +static int is_unit_row(Value *row, int pos, int len) +{ + if (!value_one_p(row[pos]) && !value_mone_p(row[pos])) + return 0; + return First_Non_Zero(row+pos+1, len-(pos+1)) == -1; +} + +/* Transform the constraints of P to "standard form". + * In particular, if P is described by + * A x + b(p) >= 0 + * then this function returns a matrix H = A U, A = H Q, such + * that D x' = D Q x >= -b(p), with D a diagonal matrix with + * positive entries. The calling function can then construct + * the standard form H' x' - I s + b(p) = 0, with H' the rows of H + * that are not positive multiples of unit vectors + * (since those correspond to D x' >= -b(p)). + * The number of rows in H' is returned in *rows_p. + * Optionally returns the matrix that maps the new variables + * back to the old variables x = U x'. + * Note that the rows of H (and P) may be reordered by this function. + */ +Matrix *standard_constraints(Polyhedron *P, unsigned nparam, int *rows_p, + Matrix **T) +{ + unsigned nvar = P->Dimension - nparam; + int i, j, d; + int rows; + Matrix *M; + Matrix *H, *U, *Q; + + assert(P->NbEq == 0); + + /* move constraints only involving parameters down + * and move unit vectors (if there are any) to the right place. + */ + for (d = 0, j = P->NbConstraints; d < j; ++d) { + int index; + index = First_Non_Zero(P->Constraint[d]+1, nvar); + if (index != -1) { + if (index != d && + is_unit_row(P->Constraint[d]+1, index, nvar)) { + Vector_Exchange(P->Constraint[d], P->Constraint[index], + P->Dimension+2); + if (index > d) + --d; + } + continue; + } + while (d < --j && First_Non_Zero(P->Constraint[j]+1, nvar) == -1) + ; + if (d >= j) + break; + Vector_Exchange(P->Constraint[d], P->Constraint[j], P->Dimension+2); + } + M = Matrix_Alloc(d, nvar); + for (j = 0; j < d; ++j) + Vector_Copy(P->Constraint[j]+1, M->p[j], nvar); + + neg_left_hermite(M, &H, &Q, &U); + Matrix_Free(M); + Matrix_Free(Q); + if (T) + *T = U; + else + Matrix_Free(U); + + /* Rearrange rows such that top of H is lower diagonal and + * compute the number of non (multiple of) unit-vector rows. + */ + rows = H->NbRows-nvar; + for (i = 0; i < H->NbColumns; ++i) { + for (j = i; j < H->NbRows; ++j) + if (value_notzero_p(H->p[j][i])) + break; + if (j != i) { + Vector_Exchange(P->Constraint[i], P->Constraint[j], P->Dimension+2); + Vector_Exchange(H->p[i], H->p[j], H->NbColumns); + } + if (First_Non_Zero(H->p[i], i) != -1) + rows++; + } + *rows_p = rows; + + return H; +} diff --git a/normalization.h b/normalization.h new file mode 100644 index 0000000..0d52c36 --- /dev/null +++ b/normalization.h @@ -0,0 +1,12 @@ +#include + +#if defined(__cplusplus) +extern "C" { +#endif + +Matrix *standard_constraints(Polyhedron *P, unsigned nparam, int *rows_p, + Matrix **T); + +#if defined(__cplusplus) +} +#endif diff --git a/polysign.c b/polysign.c index f3e83fc..c4e97b5 100644 --- a/polysign.c +++ b/polysign.c @@ -4,92 +4,6 @@ #include "polysign.h" #include "config.h" -static int is_unit_row(Value *row, int pos, int len) -{ - if (!value_one_p(row[pos]) && !value_mone_p(row[pos])) - return 0; - return First_Non_Zero(row+pos+1, len-(pos+1)) == -1; -} - -/* Transform the constraints of P to "standard form". - * In particular, if P is described by - * A x + b(p) >= 0 - * then this function returns a matrix H = A U, A = H Q, such - * that D x' = D Q x >= -b(p), with D a diagonal matrix with - * positive entries. The calling function can then construct - * the standard form H' x' - I s + b(p) = 0, with H' the rows of H - * that are not positive multiples of unit vectors - * (since those correspond to D x' >= -b(p)). - * The number of rows in H' is returned in *rows_p. - * Optionally returns the matrix that maps the new variables - * back to the old variables x = U x'. - * Note that the rows of H (and P) may be reordered by this function. - */ -Matrix *standard_constraints(Polyhedron *P, unsigned nparam, int *rows_p, - Matrix **T) -{ - unsigned nvar = P->Dimension - nparam; - int i, j, d; - int rows; - Matrix *M; - Matrix *H, *U, *Q; - - assert(P->NbEq == 0); - - /* move constraints only involving parameters down - * and move unit vectors (if there are any) to the right place. - */ - for (d = 0, j = P->NbConstraints; d < j; ++d) { - int index; - index = First_Non_Zero(P->Constraint[d]+1, nvar); - if (index != -1) { - if (index != d && - is_unit_row(P->Constraint[d]+1, index, nvar)) { - Vector_Exchange(P->Constraint[d], P->Constraint[index], - P->Dimension+2); - if (index > d) - --d; - } - continue; - } - while (d < --j && First_Non_Zero(P->Constraint[j]+1, nvar) == -1) - ; - if (d >= j) - break; - Vector_Exchange(P->Constraint[d], P->Constraint[j], P->Dimension+2); - } - M = Matrix_Alloc(d, nvar); - for (j = 0; j < d; ++j) - Vector_Copy(P->Constraint[j]+1, M->p[j], nvar); - - neg_left_hermite(M, &H, &Q, &U); - Matrix_Free(M); - Matrix_Free(Q); - if (T) - *T = U; - else - Matrix_Free(U); - - /* Rearrange rows such that top of H is lower diagonal and - * compute the number of non (multiple of) unit-vector rows. - */ - rows = H->NbRows-nvar; - for (i = 0; i < H->NbColumns; ++i) { - for (j = i; j < H->NbRows; ++j) - if (value_notzero_p(H->p[j][i])) - break; - if (j != i) { - Vector_Exchange(P->Constraint[i], P->Constraint[j], P->Dimension+2); - Vector_Exchange(H->p[i], H->p[j], H->NbColumns); - } - if (First_Non_Zero(H->p[i], i) != -1) - rows++; - } - *rows_p = rows; - - return H; -} - #ifndef HAVE_LIBGLPK enum order_sign glpk_polyhedron_affine_sign(Polyhedron *D, Matrix *T, struct barvinok_options *options) diff --git a/polysign.h b/polysign.h index 3176f7d..1d62f05 100644 --- a/polysign.h +++ b/polysign.h @@ -6,9 +6,6 @@ extern "C" { struct barvinok_options; -Matrix *standard_constraints(Polyhedron *P, unsigned nparam, int *rows_p, - Matrix **T); - enum order_sign { order_lt, order_le, order_eq, order_ge, order_gt, order_unknown, order_undefined }; diff --git a/topcom.c b/topcom.c index 53e68bb..00165f2 100644 --- a/topcom.c +++ b/topcom.c @@ -2,7 +2,7 @@ #include #include #include -#include "polysign.h" +#include "normalization.h" #include "topcom.h" #include "config.h" -- 2.11.4.GIT