1############################################################################# 2## 3#W affine-extra-4ti2gap.gi 4#W Manuel Delgado <mdelgado@fc.up.pt> 5#W Pedro Garcia-Sanchez <pedro@ugr.es> 6## 7#Y Copyright 2015-- Centro de Matemática da Universidade do Porto, Portugal and Universidad de Granada, Spain 8############################################################################# 9InstallOtherMethod(DegreesOfPrimitiveElementsOfAffineSemigroup, 10 "Computes the set of primitive elements of an affine semigroup", 11 [IsAffineSemigroup],4, 12 function(a) 13 local matrix, facs, mat, trunc, ls; 14 15 ls:=MinimalGenerators(a); 16 17 Info(InfoNumSgps,2,"Using 4ti2gap for Graver."); 18 19 mat:=TransposedMat(ls); 20 matrix := GraverBasis4ti2(["mat",mat]); 21 22 trunc:=function(ls) 23 return List(ls, y->Maximum(y,0)); 24 end; 25 26 matrix:=Set(matrix,trunc); 27 return Union(Set(matrix, x->x*ls),ls); 28end); 29 30 31InstallOtherMethod(HilbertBasisOfSystemOfHomogeneousEquations, 32 "Computes a Hilbert basis of a systemd of linear Diophantine equations, some eventually in congruences.", 33 [IsRectangularTable,IsHomogeneousList],6, 34 function(ls,md) 35 local homogeneous, withCongruences; 36 37 homogeneous:= function(l) 38 local problem, matrix,mat,sign; 39 40 Info(InfoNumSgps,2,"Using 4ti2gap for Hilbert."); 41 42 #if not(IsRectangularTable(l)) then 43 # Error("The argument must be a matrix."); 44 #fi; 45 if not(IsInt(l[1][1])) then 46 Error("The matrix must be of integers."); 47 fi; 48 49 mat:=l; 50 sign:=[List(l[1],_->1)]; 51 problem:=["mat",mat, "sign", sign]; 52 53 matrix := HilbertBasis4ti2(problem).zhom; 54 return matrix; 55 56 end; 57 58 withCongruences:=function(ls,md) 59 local l,n,m,diag,dim,d, hil, zero, leq; 60 61 leq:= function(v1,v2) 62 local v; 63 v:=v2-v1; 64 return (First(v,n->n<0)=fail); 65 end; 66 67 #if not(IsRectangularTable(ls)) then 68 # Error("The first argument must be a matrix."); 69 #fi; 70 71 if not(IsListOfIntegersNS(md)) or ForAny(md, x->not(IsPosInt(x))) then 72 Error("The second argument must be a list of positive integers."); 73 fi; 74 75 n:=Length(ls); 76 dim:=Length(ls[1]); 77 m:=Length(md); 78 if m>n then 79 Error("There are more modulus than equations."); 80 fi; 81 82 diag:=Concatenation(md,List([1..n-m],_->0)); 83 d:=DiagonalMat(diag); 84 l:=TransposedMat(Concatenation(TransposedMat(ls),d,-d)); 85 zero:=List([1..dim],_->0); 86 87 hil:=Difference(List(homogeneous(l), x->x{[1..dim]}),[zero]); 88 return hil; 89 90 return Filtered(hil, y->Filtered(hil,x->leq(x,y))=[y]); 91 end; 92 ## end of local functions ... 93 94 #ls := arg[1][1]; 95 #md := arg[1][2]; 96 if md = [] then 97 return homogeneous(ls); 98 else 99 return withCongruences(ls,md); 100 101 fi; 102 103end); 104 105InstallOtherMethod(HilbertBasisOfSystemOfHomogeneousInequalities, 106 "Computes a Hilbert basis of l*x>=0, x>=0", 107 [IsRectangularTable],6, 108 function(l) 109 local problem, matrix,mat,sign,rel; 110 111 Info(InfoNumSgps,2,"Using 4ti2gap for Hilbert."); 112 113 #if not(IsRectangularTable(l)) then 114 # Error("The argument must be a matrix."); 115 #fi; 116 if not(IsInt(l[1][1])) then 117 Error("The matrix must be of integers."); 118 fi; 119 120 mat:=l; 121 sign:=[List(l[1],_->1)]; 122 rel:=[List(l[1],_->">")]; 123 problem:=["mat",mat,"rel",rel,"sign",sign]; 124 matrix:=HilbertBasis4ti2(problem); 125 return matrix; 126 127end); 128 129 130InstallOtherMethod(FactorizationsVectorWRTList, 131 "Computes the factorizations of v in terms of the elments in ls", 132 [IsHomogeneousList,IsMatrix],6, 133 function(v,l) 134 local matrix,mat,rhs,sign,problem, n; 135 136 Info(InfoNumSgps,2,"Using 4ti2gap for factorizations."); 137 138 if not(IsListOfIntegersNS(v)) then 139 Error("The first argument must be a list of integers."); 140 fi; 141 142 if not(IsInt(l[1][1])) then 143 Error("The matrix must be of integers."); 144 fi; 145 146 mat:=TransposedMat(Concatenation(l,[-v])); 147 if not(IsRectangularTable(mat)) then 148 Error("The list in the second argument must have the same length as all the lists in the first argument."); 149 fi; 150 151 sign:=[List(l,_->1)]; 152 rhs:=[v]; 153 problem:=["mat",TransposedMat(l),"sign",sign,"rhs",rhs]; 154 matrix := ZSolve4ti2(problem); 155 return matrix.zinhom; 156 157end); 158 159 160InstallOtherMethod(GeneratorsOfKernelCongruence, 161 "Computes a set of generators of the kernel congruence of the monoid morphism associated to a matrix", 162 [IsRectangularTable],7, 163 function(m) 164 local positivenegative, gr; 165 166 positivenegative:=function(p) 167 local d1, d2; 168 d1:=List(p, i->Maximum(i,0)); 169 d2:=List(p, i->-Minimum(0,i)); 170 return [d1,d2]; 171 end; 172 173 if not(ForAll(m, l->ForAll(l, x->(x=0) or IsPosInt(x)))) then 174 Error("The argument must be a matrix of nonnegative integer."); 175 fi; 176 177 gr:=GroebnerBasis4ti2(TransposedMat(m)); 178 Info(InfoNumSgps,2,"4ti output:",gr); 179 180 return List(gr, x->positivenegative(x)); 181end); 182 183############################################################ 184# computes a canonical basis of the kernel congruence 185# of the monoid morphism associated to the matrix m with 186# nonnegative integer coefficients wrt the term ordering 187# the kernel is the pairs (x,y) such that xm=ym 188############################################################ 189InstallMethod(CanonicalBasisOfKernelCongruence, 190"Computes a canonical basis for the congruence of of the monoid morphism associated to the matrix", 191 [IsRectangularTable, IsMonomialOrdering],7, 192 function(m,ord) 193 local positivenegative, gr, nord, to; 194 195 positivenegative:=function(p) 196 local d1, d2; 197 d1:=List(p, i->Maximum(i,0)); 198 d2:=List(p, i->-Minimum(0,i)); 199 return [d1,d2]; 200 end; 201 202 if not(ForAll(m, l->ForAll(l, x->(x=0) or IsPosInt(x)))) then 203 Error("The argument must be a matrix of nonnegative integer."); 204 fi; 205 206 207 # trick taken from the package Singular 208 nord := Name( ord ); 209 nord := nord{[ 1 .. Position( nord, '(' ) - 1 ]}; 210 if nord = "MonomialLexOrdering" then 211 to := "lex"; 212 elif nord = "MonomialGrevlexOrdering" then 213 to := "grevlex"; 214 elif nord = "MonomialGrlexOrdering" then 215 to := "grlex"; 216 else 217 Error( "the ordering ", ord, " is not yet supported\n" ); 218 fi; 219 220 gr:=GroebnerBasis4ti2(TransposedMat(m),to); 221 Info(InfoNumSgps,2,"4ti output:",gr); 222 223 return Set(gr, x->positivenegative(x)); 224 end); 225 226############################################################ 227# computes the Graver basis of matrix with integer entries 228############################################################ 229InstallMethod(GraverBasis, 230 "Computes the Graver basis of the matrix", 231 [IsRectangularTable],8, 232 function(a) 233 #4ti2gap implementation 234 local gr; 235 236 237 if not(IsRectangularTable(a)) then 238 Error("The argument must be a matrix."); 239 fi; 240 241 if not(IsInt(a[1][1])) then 242 Error("The entries of the matrix must be integers."); 243 fi; 244 245 Info(InfoNumSgps,2,"Using 4ti2gap for Graver."); 246 247 248 gr := GraverBasis4ti2(["mat",a]); 249 return Union(gr,-gr); 250 end); 251 252 253 254InstallOtherMethod(MinimalPresentationOfAffineSemigroup, 255 "Computes a minimimal presentation of the affine semigroup", 256 [IsAffineSemigroup],6, 257 function(a) 258 local gens, positive, gr, candidates, pres, rclass,exps, c; 259 260 positive:=function(x) 261 local p,i; 262 263 p:=[]; 264 for i in [1..Length(x)] do 265 p[i]:=Maximum(x[i],0); 266 od; 267 268 return p; 269 end; 270 if not(IsAffineSemigroup(a)) then 271 Error("The argument must be an affine semigroup."); 272 fi; 273 274 gens:=MinimalGenerators(a); 275 276 gr:=GroebnerBasis4ti2(TransposedMat(gens)); 277 Info(InfoNumSgps,2,"4ti output:",gr); 278 279 candidates:=Set(gr,q->positive(q)); 280 candidates:=Set(candidates,c->c*gens); 281 Info(InfoNumSgps,2, "Candidates to Betti elements",candidates); 282 pres:=[]; 283 for c in candidates do 284 exps:=FactorizationsVectorWRTList(c,gens); 285 rclass:=RClassesOfSetOfFactorizations(exps); 286 if Length(rclass)>1 then 287 pres:=Concatenation(pres,List([2..Length(rclass)], 288 i->[rclass[1][1],rclass[i][1]])); 289 fi; 290 od; 291 return pres; 292end); 293 294 295 296##################################################################### 297# Computes the omega-primality of v in the affine semigroup a 298##################################################################### 299InstallOtherMethod(OmegaPrimalityOfElementInAffineSemigroup, 300 "Computes the omega-primality of v in the affine semigroup a", 301 [IsHomogeneousList,IsAffineSemigroup],6, 302 function(v,a) 303 local ls, n, mat,extfact,par,tot,le; 304 305 le:=function(a,b) #ordinary partial order 306 return ForAll(b-a,x-> x>=0); 307 end; 308 309 if not(IsAffineSemigroup(a)) then 310 Error("The second argument must be an affine semigroup"); 311 fi; 312 313 if not(IsListOfIntegersNS(v)) then 314 Error("The first argument must be a list of integers."); 315 fi; 316 317 if not(ForAll(v, x-> x>=0)) then 318 Error("The first argument must be a list of on nonnegative integers."); 319 fi; 320 321 ls:=MinimalGenerators(a); 322 n:=Length(ls); 323 mat:=TransposedMat(Concatenation(ls,-ls,[-v])); 324 325 if not(IsRectangularTable(mat)) then 326 Error("The first argument has not the dimension of the second."); 327 fi; 328 329 extfact:=FactorizationsVectorWRTList(v,Concatenation(ls,-ls)); 330 331 par:=Set(extfact, f->f{[1..n]}); 332 tot:=Filtered(par, f-> Filtered(par, g-> le(g,f))=[f]); 333 Info(InfoNumSgps,2,"Minimals of v+ls =",tot); 334 if tot=[] then 335 return 0; 336 fi; 337 return Maximum(Set(tot, Sum)); 338end); 339 340##################################################################### 341# Computes the omega-primality of v in the full affine semigroup a 342##################################################################### 343InstallOtherMethod(OmegaPrimalityOfElementInAffineSemigroup, 344 "Computes the omega-primality of v in the full affine semigroup a", 345 [IsHomogeneousList,IsAffineSemigroup and HasEquations],6, 346 function(v,a) 347 local ls, n, mat,extfact,par,tot,le; 348 349 le:=function(a,b) #ordinary partial order 350 return ForAll(b-a,x-> x>=0); 351 end; 352 353 if not(IsAffineSemigroup(a)) then 354 Error("The second argument must be an affine semigroup"); 355 fi; 356 357 if not(IsListOfIntegersNS(v)) then 358 Error("The first argument must be a list of integers."); 359 fi; 360 361 if not(ForAll(v, x-> x>=0)) then 362 Error("The first argument must be a list of on nonnegative integers."); 363 fi; 364 365 ls:=MinimalGenerators(a); 366 n:=Length(ls); 367 mat:=TransposedMat(Concatenation(ls,-ls,[-v])); 368 369 if not(IsRectangularTable(mat)) then 370 Error("The first argument has not the dimension of the second."); 371 fi; 372 373 Info(InfoNumSgps,2,"Using 4ti2gap with full affine semigroup"); 374 375 extfact:=ZSolve4ti2(["mat",TransposedMat(ls),"rel",List(v,_->1), 376 "sign",List([1..n],_->1),"rhs",v ]); 377 378 tot:=extfact.zinhom; 379 Info(InfoNumSgps,2,"Minimals of v+ls =",tot); 380 if tot=[] then 381 return 0; 382 fi; 383 return Maximum(Set(tot, Sum)); 384end); 385 386#ZSolve4ti2(["mat",TransposedMat([[2,0],[0,2],[1,2],[2,1]]),"rel",[1,1],"sign",[1,1,1,1],"rhs",[[15,15]]]); 387 388######## 389# Tame degree for full affine semigroups 390######## 391InstallMethod(TameDegreeOfAffineSemigroup, 392 "Computes the tame degree of the full affine semigroup a", 393 [IsAffineSemigroup and HasEquations],2, 394 function(a) 395 local ls, min, tame, gen,m,n, facts, t, minfacts; 396 397 Info(InfoNumSgps,2,"Using 4ti2gap with full affine semigroup"); 398 399 ls:=MinimalGenerators(a); 400 tame:=0; 401 n:=Length(ls); 402 403 for gen in ls do 404 minfacts:=ZSolve4ti2(["mat",TransposedMat(ls),"rel",List(gen,_->1), 405 "sign",List([1..n],_->1),"rhs",gen ]).zinhom; 406 min:=List(minfacts, x->x*ls); 407 Info(InfoNumSgps,2,"Minimal elements of ",gen,"+a=",min); 408 for m in min do 409 facts:=FactorizationsVectorWRTList(m,ls); 410 t:=TameDegreeOfSetOfFactorizations(facts); 411 if t> tame then 412 tame:=t; 413 Info(InfoNumSgps,2,"Tame degree updated to ",tame); 414 fi; 415 od; 416 od; 417 return tame; 418 419end); 420