diff --git a/docs/oscar_references.bib b/docs/oscar_references.bib index e7963429c25b..ec4cd965b01f 100644 --- a/docs/oscar_references.bib +++ b/docs/oscar_references.bib @@ -37,6 +37,30 @@ @Article{AG10 eprint = {0810.1148} } +@InProceedings{AGK96, + author = {Amrhein, Beatrice and Gloor, Oliver and Küchlin, Wolfgang}, + title = {Walking faster}, + booktitle = {Design and Implementation of Symbolic Computation Systems}, + series = {DISCO 1996}, + address = {Heidelberg, Germany}, + publisher = {Springer Berlin, Heidelberg}, + pages = {150--161}, + year = {1996}, + doi = {10.1007/3-540-61697-7_14}, + location = {Karlsruhe, Germany} +} + +@Article{AGK97, + author = {Amrhein, Beatrice and Gloor, Oliver and Küchlin, Wolfgang}, + title = {On the walk}, + journal = {Theoretical Computer Science}, + volume = {187}, + number = {1}, + pages = {179--202}, + year = {1997}, + doi = {10.1016/s0304-3975(97)00064-9} +} + @Article{AHK18, author = {Adiprasito, Karim and Huh, June and Katz, Eric}, title = {Hodge theory for combinatorial geometries}, @@ -453,6 +477,27 @@ @Article{CHM98 doi = {10.1112/S1461157000000115} } +@Article{CKM97, + author = {Collart, S. and Kalkbrener, M. and Mall, D.}, + title = {Converting Bases with the Gröbner Walk}, + journal = {Journal of Symbolic Computation}, + volume = {24}, + number = {3}, + pages = {465--469}, + year = {1997}, + doi = {10.1006/jsco.1996.0145} +} + +@Book{CLO05, + author = {Cox, David A. and Little, John and O'Shea, Donal}, + title = {Using Algebraic Geometry}, + series = {Graduate Texts in Mathematics}, + volume = {185}, + publisher = {Springer-Verlag}, + year = {2005}, + doi = {10.1007/b138611} +} + @Book{CLS11, author = {Cox, David A. and Little, John B. and Schenck, Henry K.}, title = {Toric varieties}, @@ -910,6 +955,17 @@ @Article{FGLM93 doi = {10.1006/jsco.1993.1051} } +@Article{FJLT07, + author = {Fukuda, K. and Jensen, A. N. and Lauritzen, N. and Thomas, R.}, + title = {The generic Gröbner walk}, + journal = {Journal of Symbolic Computation}, + volume = {42}, + number = {3}, + pages = {298--312}, + year = {2007}, + doi = {10.1016/j.jsc.2006.09.004} +} + @Article{FJR17, author = {Fan, Huijun and Jarvis, Tyler and Ruan, Yongbin}, title = {A mathematical theory of the gauged linear sigma model}, @@ -923,6 +979,18 @@ @Article{FJR17 doi = {10.2140/gt.2018.22.235} } +@Article{FJT07, + author = {Fukuda, Komei and Jensen, Anders N. and Thomas, Rekha R.}, + title = {Computing Groebner Fans}, + journal = {Math. Comp.}, + fjournal = {Mathematics of Computation}, + volume = {76}, + number = {260}, + pages = {2189--2213}, + year = {2007}, + doi = {10.1090/S0025-5718-07-01986-2} +} + @InProceedings{FLINT, bibkey = {FLINT}, author = {W. B. Hart}, @@ -2161,6 +2229,28 @@ @Article{Tay87 primaryclass = {math.GR} } +@Article{Tra00, + author = {Tran, Quoc-Nam}, + title = {A Fast Algorithm for Gröbner Basis Conversion and its Applications}, + journal = {Journal of Symbolic Computation}, + volume = {30}, + number = {4}, + pages = {451--467}, + year = {2000}, + doi = {10.1006/jsco.1999.0416} +} + +@Article{Tra04, + author = {Tran, Quoc-Nam}, + title = {Efficient Groebner walk conversion for implicitization of geometric objects}, + journal = {Computer Aided Geometric Design}, + volume = {21}, + number = {9}, + pages = {837--857}, + year = {2004}, + doi = {10.1016/j.cagd.2004.07.001} +} + @Misc{VE22, author = {Vaughan-Lee, Michael and Eick, Bettina}, title = {SglPPow, Database of groups of prime-power order for some prime-powers, Version 2.3}, diff --git a/experimental/GroebnerWalk/README.md b/experimental/GroebnerWalk/README.md new file mode 100644 index 000000000000..b9132345d236 --- /dev/null +++ b/experimental/GroebnerWalk/README.md @@ -0,0 +1,49 @@ +# Gröbner walk + +[![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.11065978.svg)](https://doi.org/10.5281/zenodo.11065978) + +`GroebnerWalk` provides implementations of Gröbner walk algorithms +for computing Gröbner bases over fields on top of Oscar.jl. + +## Usage + +This module provides the function `groebner_walk` as interface to the algorithms. +The following example demonstrates the usage. First, we define the ideal Oscar.jl. +```julia +using Oscar + +R, (x,y) = QQ[:x, :y] # define ring ... +I = ideal([y^4+ x^3-x^2+x,x^4]) # ... and ideal +``` +Then, we can pass the ideal to `groebner_walk` to calculate the Gröbner basis. + +By default, `groebner_walk` starts with a Gröbner basis with respect to the default ordering on `R` +and converts this into a Gröbner basis with respect to the lexicographic ordering on `R`. +This is what the following code block accomplishes. +```julia +using Oscar + +groebner_walk(I) # compute the Groebner basis +``` +If one wants to specify `target` and `start` orderings explicitly, above function call needs to be written as follows. +```julia +groebner_walk(I, lex(R), default_ordering(R)) # compute the Groebner basis +``` + +Additionally, there are certain special ideals provided that are used for benchmarking +of this module. + +## Status +At the moment, the standard walk by Collart, Kalkbrener and Mall (1997) and the generic walk by Fukuda et al. (2007) are implemented. + +## Contacts +The library is maintained by Kamillo Ferry (kafe (at) kafe (dot) dev) and Francesco Nowell (francesconowell (at) gmail (dot) com). + +## Acknowledgement +The current implementation is based on an implementation by Jordi Welp. We thank him for +laying the groundwork for this package. + +## References +- Collart, S., M. Kalkbrener, and D. Mall. ‘Converting Bases with the Gröbner Walk’. Journal of Symbolic Computation 24, no. 3–4 (September 1997): 465–69. https://doi.org/10.1006/jsco.1996.0145. +- Fukuda, K., A. N. Jensen, N. Lauritzen, and R. Thomas. ‘The Generic Gröbner Walk’. Journal of Symbolic Computation 42, no. 3 (1 March 2007): 298–312. https://doi.org/10.1016/j.jsc.2006.09.004. + diff --git a/experimental/GroebnerWalk/benchmark/agk.jl b/experimental/GroebnerWalk/benchmark/agk.jl new file mode 100644 index 000000000000..61d9c2a5dc99 --- /dev/null +++ b/experimental/GroebnerWalk/benchmark/agk.jl @@ -0,0 +1,17 @@ +using Oscar + +function agk4(k::Field) + R, (x,y,z,u,v) = polynomial_ring(k, ["x","y","z","u","v"]) + + o1 = weight_ordering([1,1,1,0,0], degrevlex(R)) + o2 = weight_ordering([0,0,0,1,1], degrevlex(R)) + + F = [ + u + u^2 - 2*v - 2*u^2*v + 2*u*v^2 - x, + -6*u + 2*v + v^2 - 5*v^3 + 2*u*v^2 - 4*u^2*v^2 - y, + -2 + 2*u^2 + 6*v - 3*u^2*v^2 - z + ] + + return ideal(F), o2, o1 +end + diff --git a/experimental/GroebnerWalk/benchmark/benchmarks.jl b/experimental/GroebnerWalk/benchmark/benchmarks.jl new file mode 100644 index 000000000000..6feb716153b4 --- /dev/null +++ b/experimental/GroebnerWalk/benchmark/benchmarks.jl @@ -0,0 +1,65 @@ +using Oscar + +include("cyclic.jl") +include("katsura.jl") +include("agk.jl") +include("tran3.3.jl") + +function benchmark( + io, + name::String, + I::Ideal, + target::MonomialOrdering, + start::MonomialOrdering +) + print(io, name, ","); flush(io) + t = @elapsed groebner_walk($I, $target, $start; algorithm=:standard) + print(io, t, ","); flush(io) + t = @elapsed groebner_walk($I, $target, $start; algorithm=:generic) + print(io, t, ","); flush(io) + t = @elapsed groebner_walk($I, $target, $start; algorithm=:perturbed) + print(io, t, ","); flush(io) + t = @elapsed groebner_basis($I; ordering=$target) + println(io, t); flush(io) +end + +function print_header(io) + print(io, "name,standard_walk,generic_walk,perturbed_walk,buchberger\n") +end + +p = 11863279 +Fp = GF(p) +open("results.csv", "a") do io + R, (x, y) = QQ[:x, :y] + I = ideal([y^4 + x^3 - x^2 + x, x^4]) + benchmark(io, "simple", I, lex(R), default_ordering(R)) + + benchmark(io, "cyclic5-QQ", cyclic5(QQ)...) + benchmark(io, "cyclic5-Fp", cyclic5(Fp)...) + + benchmark(io, "cyclic6-QQ", cyclic6(QQ)...) + benchmark(io, "cyclic6-Fp", cyclic6(Fp)...) + + I = katsura(6) + R = base_ring(I) + benchmark(io, "katsura6-QQ", I, lex(R), default_ordering(R)) + + Ip = map(gens(I)) do f + change_coefficient_ring(Fp, f) + end |> ideal + Rp = base_ring(Ip) + benchmark(io, "katsura6-Fp", Ip, lex(Rp), default_ordering(Rp)) + + benchmark(io, "cyclic7-QQ", cyclic7(QQ)...) + benchmark(io, "cyclic7-Fp", cyclic7(Fp)...) + + benchmark(io, "agk4-QQ", agk4(QQ)...) + benchmark(io, "agk4-Fp", agk4(Fp)...) + + benchmark(io, "tran3.3-QQ", tran33(QQ)...) + benchmark(io, "tran3.3-Fp", tran33(Fp)...) + + benchmark(io, "newellp1-QQ", Oscar.newell_patch_with_orderings(QQ)...) + benchmark(io, "newellp1-Fp", Oscar.newell_patch_with_orderings(Fp)...) +end + diff --git a/experimental/GroebnerWalk/benchmark/cyclic.jl b/experimental/GroebnerWalk/benchmark/cyclic.jl new file mode 100644 index 000000000000..4e1c829459e7 --- /dev/null +++ b/experimental/GroebnerWalk/benchmark/cyclic.jl @@ -0,0 +1,79 @@ +using Oscar + +function cyclic5(k::Field=QQ) + R, (a,b,c,d,x) = k[:a,:b,:c,:d,:x] + F = [ + a + b + c + d + x, + a*b + b*c + c*d + d*x + x*a, + a*b*c + b*c*d + c*d*x + d*x*a + x*a*b, + a*b*c*d + b*c*d*x + c*d*x*a + d*x*a*b + x*a*b*c, + a*b*c*d*x - 1 + ] + + return ideal(F), lex(R), default_ordering(R) +end + +function cyclic6(k::Field=QQ) + R, (z0,z1,z2,z3,z4,z5) = k[:z0,:z1,:z2,:z3,:z4,:z5] + F = [ + z0 + z1 + z2 + z3 + z4 + z5, + z0*z1 + z1*z2 + z2*z3 + z3*z4 + z4*z5 + z5*z0, + z0*z1*z2 + z1*z2*z3 + z2*z3*z4 + z3*z4*z5 + z4*z5*z0 + z5*z0*z1, + z0*z1*z2*z3 + z1*z2*z3*z4 + z2*z3*z4*z5 + z3*z4*z5*z0 + z4*z5*z0*z1 + + z5*z0*z1*z2, + z0*z1*z2*z3*z4 + z1*z2*z3*z4*z5 + z2*z3*z4*z5*z0 + z3*z4*z5*z0*z1 + + z4*z5*z0*z1*z2 + z5*z0*z1*z2*z3, + z0*z1*z2*z3*z4*z5 - 1 + ] + + return ideal(F), lex(R), default_ordering(R) +end + +function cyclic7(k::Field=QQ) + R, (z0, z1, z2, z3, z4, z5, z6) = QQ[:z0, :z1, :z2, :z3, :z4, :z5, :z6] + F = [ + z0 + z1 + z2 + z3 + z4 + z5 + z6, + z0 * z1 + z1 * z2 + z2 * z3 + z3 * z4 + z4 * z5 + z5 * z6 + z6 * z0, + z0 * z1 * z2 + z1 * z2 * z3 + z2 * z3 * z4 + z3 * z4 * z5 + z4 * z5 * z6 + z5 * z6 * z0 + z6 * z0 * z1, + z0 * z1 * z2 * z3 + z1 * z2 * z3 * z4 + z2 * z3 * z4 * z5 + z3 * z4 * z5 * z6 + z4 * z5 * z6 * z0 + + z5 * z6 * z0 * z1 + z6 * z0 * z1 * z2, + z0 * z1 * z2 * z3 * z4 + z1 * z2 * z3 * z4 * z5 + z2 * z3 * z4 * z5 * z6 + z3 * z4 * z5 * z6 * z0 + + z4 * z5 * z6 * z0 * z1 + z5 * z6 * z0 * z1 * z2 + z6 * z0 * z1 * z2 * z3, + z0 * z1 * z2 * z3 * z4 * z5 + z1 * z2 * z3 * z4 * z5 * z6 + z2 * z3 * z4 * z5 * z6 * z0 + z3 * z4 * z5 * z6 * z0 * z1 + + z4 * z5 * z6 * z0 * z1 * z2 + z5 * z6 * z0 * z1 * z2 * z3 + z6 * z0 * z1 * z2 * z3 * z4, + z0 * z1 * z2 * z3 * z4 * z5 * z6 - 1 + ] + + return ideal(F), lex(R), default_ordering(R) +end + +function cyclic8(k::Field=QQ) + R, (z0, z1, z2, z3, z4, z5, z6, z7) = k[:z0, :z1, :z2, :z3, :z4, :z5, :z6, :z7] + + F = [ + z0 + z1 + z2 + z3 + z4 + z5 + z6 + z7, + + z0*z1 + z1*z2 + z2*z3 + z3*z4 + z4*z5 + z5*z6 + z6*z7 + z7*z0, + + z0*z1*z2 + z1*z2*z3 + z2*z3*z4 + z3*z4*z5 + z4*z5*z6 + z5*z6*z7 + + z6*z7*z0 + z7*z0*z1, + + z0*z1*z2*z3 + z1*z2*z3*z4 + z2*z3*z4*z5 + z3*z4*z5*z6 + z4*z5*z6*z7 + + z5*z6*z7*z0 + z6*z7*z0*z1 + z7*z0*z1*z2, + + z0*z1*z2*z3*z4 + z1*z2*z3*z4*z5 + z2*z3*z4*z5*z6 + z3*z4*z5*z6*z7 + + z4*z5*z6*z7*z0 + z5*z6*z7*z0*z1 + z6*z7*z0*z1*z2 + z7*z0*z1*z2*z3, + + z0*z1*z2*z3*z4*z5 + z1*z2*z3*z4*z5*z6 + z2*z3*z4*z5*z6*z7 + z3*z4*z5*z6*z7*z0 + + z4*z5*z6*z7*z0*z1 + z5*z6*z7*z0*z1*z2 + z6*z7*z0*z1*z2*z3 + z7*z0*z1*z2*z3*z4, + + z0*z1*z2*z3*z4*z5*z6 + z1*z2*z3*z4*z5*z6*z7 + z2*z3*z4*z5*z6*z7*z0 + + z3*z4*z5*z6*z7*z0*z1 + z4*z5*z6*z7*z0*z1*z2 + z5*z6*z7*z0*z1*z2*z3 + + z6*z7*z0*z1*z2*z3*z4 + z7*z0*z1*z2*z3*z4*z5, + + z0*z1*z2*z3*z4*z5*z6*z7 - 1 + ] + + return ideal(F), lex(R), default_ordering(R) +end + diff --git a/experimental/GroebnerWalk/benchmark/tran3.3.jl b/experimental/GroebnerWalk/benchmark/tran3.3.jl new file mode 100644 index 000000000000..565a108a482a --- /dev/null +++ b/experimental/GroebnerWalk/benchmark/tran3.3.jl @@ -0,0 +1,10 @@ +using Oscar + +function tran33(k::Field) + R, (x,y,z) = k[:x,:y,:z] + + F = [16 + 3*x^3+16*x^2*z+14*x^2*y^3, 6+y^3*z+17*x^2*z^2+7*x*y^2*z^2+13*x^3*z^2] + + return ideal(F), lex(R), default_ordering(R) +end + diff --git a/experimental/GroebnerWalk/docs/doc.main b/experimental/GroebnerWalk/docs/doc.main new file mode 100644 index 000000000000..0264d8306ed9 --- /dev/null +++ b/experimental/GroebnerWalk/docs/doc.main @@ -0,0 +1,6 @@ +[ + "Gröbner walk" => [ + "introduction.md", + "special-ideals.md" + ] +] diff --git a/experimental/GroebnerWalk/docs/src/introduction.md b/experimental/GroebnerWalk/docs/src/introduction.md new file mode 100644 index 000000000000..2583cb9056cd --- /dev/null +++ b/experimental/GroebnerWalk/docs/src/introduction.md @@ -0,0 +1,22 @@ +```@meta +CurrentModule = Oscar +``` + +# Usage + +The Gröbner walk is an approach to reduce the computational complexity of Gröbner basis computations as proposed by [AGK97](@cite). +These incarnations of the Gröbner walk refer to a family of algorithms that perform a reverse local search on the cones of the Gröbner fan. +Then, a Gröbner basis is calculated for each encountered cone while reusing the generators obtained from the previous cone. + +The implemented algorithms may be accessed using the following function. + +```@docs + groebner_walk( + I::MPolyIdeal, + target::MonomialOrdering = lex(base_ring(I)), + start::MonomialOrdering = default_ordering(base_ring(I)); + perturbation_degree = ngens(base_ring(I)), + algorithm::Symbol = :standard + ) +``` + diff --git a/experimental/GroebnerWalk/docs/src/special-ideals.md b/experimental/GroebnerWalk/docs/src/special-ideals.md new file mode 100644 index 000000000000..6c58ef1fe75a --- /dev/null +++ b/experimental/GroebnerWalk/docs/src/special-ideals.md @@ -0,0 +1,16 @@ +```@meta +CurrentModule = Oscar +``` + +# Special ideals used for benchmarking + +We bundle a couple of special ideals useful for benchmarking of the Gröbner walk. + +```@docs + newell_patch(k::Union{QQField, QQBarFieldElem}, n::Int=1) + newell_patch(k::Field, n::Int=1) +``` + +```@docs + newell_patch_with_orderings(k::Field, n::Int=1) +``` diff --git a/experimental/GroebnerWalk/examples/implicitization/agk4.jl b/experimental/GroebnerWalk/examples/implicitization/agk4.jl new file mode 100644 index 000000000000..cfc243a63788 --- /dev/null +++ b/experimental/GroebnerWalk/examples/implicitization/agk4.jl @@ -0,0 +1,19 @@ +#= + Implicitization of Bezier surface. Example taken from Amrhein, Gloor, Küchlin. "On the walk" (2004) +=# +using Oscar + +R, (x,y,z,u,v) = polynomial_ring(QQ, ["x","y","z","u","v"]) + +o1 = matrix_ordering(R, [1 1 1 0 0; 0 0 0 1 1; 0 0 0 1 0; 1 1 0 0 0; 1 0 0 0 0]) +o2 = matrix_ordering(R, [0 0 0 1 1; 1 1 1 0 0; 1 1 0 0 0; 1 0 0 0 0; 0 0 0 1 0]) + +I = ideal([ + u + u^2 - 2*v - 2*u^2*v + 2*u*v^2 - x, + -6*u + 2*v + v^2 - 5*v^3 + 2*u*v^2 - 4*u^2*v^2 - y, + -2 + 2*u^2 + 6*v - 3*u^2*v^2 - z +]) +set_verbosity_level(:groebner_walk, 1) + +G = groebner_walk(I, o2, o1; algorithm=:standard) + diff --git a/experimental/GroebnerWalk/examples/implicitization/newellp1.jl b/experimental/GroebnerWalk/examples/implicitization/newellp1.jl new file mode 100644 index 000000000000..feaf316bf6dc --- /dev/null +++ b/experimental/GroebnerWalk/examples/implicitization/newellp1.jl @@ -0,0 +1,11 @@ +#= + Newell's teapot, patch 1 + A 2-dimensional ideal from Tran. "Efficient Groebner walk conversion for implicitization of geometric objects" (2004) +=# +using Oscar + +I, o1, o2 = Oscar.newell_patch_with_orderings(QQ, 1) + +set_verbosity_level(:groebner_walk, 1) +G = groebner_walk(I, o2, o1; algorithm=:standard) + diff --git a/experimental/GroebnerWalk/examples/implicitization/tran3.3.jl b/experimental/GroebnerWalk/examples/implicitization/tran3.3.jl new file mode 100644 index 000000000000..f5dab55b5468 --- /dev/null +++ b/experimental/GroebnerWalk/examples/implicitization/tran3.3.jl @@ -0,0 +1,11 @@ +#= + Example 3.3 from Tran. "A fast algorithm for Gröbner basis conversion and its applications" (2000) + A 1-dimensional ideal (application not given) +=# +using Oscar +R, (x,y,z) = polynomial_ring(QQ, ["x","y","z"]) +I = ideal([16 + 3*x^3+16*x^2*z+14*x^2*y^3, 6+y^3*z+17*x^2*z^2+7*x*y^2*z^2+13*x^3*z^2]) + +set_verbosity_level(:groebner_walk, 1) +Gs = groebner_walk(I, lex(R), algorithm =:standard) + diff --git a/experimental/GroebnerWalk/examples/knapsack/fukuda_cuww1.jl b/experimental/GroebnerWalk/examples/knapsack/fukuda_cuww1.jl new file mode 100644 index 000000000000..c36928d90d77 --- /dev/null +++ b/experimental/GroebnerWalk/examples/knapsack/fukuda_cuww1.jl @@ -0,0 +1,23 @@ +#= + Hard integer knapsack problem from Aardal, Karen and Lenstra. ‘Hard Equality Constrained Integer Knapsacks’. (2004) +=# + +using Oscar +R, (t, x1, x2, x3, x4, x5) = polynomial_ring(QQ, ["t", "x1", "x2", "x3", "x4", "x5"]) + +f = 12223 * x1 + 12224 * x2 + 36674 * x3 + 61119 * x4 + 85569 * x5 - 89643481 + +f1 = x1 - t^1223 +f2 = x2 - t^1224 +f3 = x3 - t^36674 +f4 = x4 - t^61119 +f5 = x5 - t^85569 +I = ideal([f1, f2, f3, f4, f5]) + +set_verbosity_level(:groebner_walk, 1) + +o_t = weight_ordering([1, 0, 0, 0, 0, 0], degrevlex(R)) +o_s = weight_ordering([0, 1, 1, 1, 1, 1], degrevlex(R)) + +ts = @elapsed Gs = groebner_walk(I, o_t, o_s, algorithm=:standard) + diff --git a/experimental/GroebnerWalk/examples/knapsack/fukuda_prob6.jl b/experimental/GroebnerWalk/examples/knapsack/fukuda_prob6.jl new file mode 100644 index 000000000000..fe65605a33f6 --- /dev/null +++ b/experimental/GroebnerWalk/examples/knapsack/fukuda_prob6.jl @@ -0,0 +1,45 @@ +using Oscar +R, (t, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10) = polynomial_ring(QQ, ["t","x1","x2", "x3", "x4", "x5", "x6", "x7", "x8", "x9", "x10" ]) + +f = 12223*x1 + 12224*x2 +36674*x3+61119*x4+85569*x5 -89643481 + +f1 = x1 - t^20601 +f2 = x2 - t^40429 +f3 = x3 - t^42207 +f4 = x4 - t^45415 +f5 = x5 - t^53725 +f6 = x6 - t^61919 +f7 = x7 - t^64670 +f8 = x8 - t^69340 +f9 = x9 - t^78539 +f10 = x10 - t^95043 +I = ideal([f1, f2, f3, f4, f5, f6, f7, f8, f9, f10]) + +set_verbosity_level(:groebner_walk, 1) +o_t = matrix_ordering(R, [ +1 0 0 0 0 0 0 0 0 0 0; +0 1 1 1 1 1 1 1 1 1 1; +0 1 1 1 1 1 1 1 1 1 0; +0 1 1 1 1 1 1 1 1 0 0; +0 1 1 1 1 1 1 1 0 0 0; +0 1 1 1 1 1 1 0 0 0 0; +0 1 1 1 1 1 0 0 0 0 0; +0 1 1 1 1 0 0 0 0 0 0; +0 1 1 1 0 0 0 0 0 0 0; +0 1 1 0 0 0 0 0 0 0 0; +0 1 0 0 0 0 0 0 0 0 0]) +o_s = matrix_ordering(R, [ +0 1 1 1 1 1 1 1 1 1 1; +1 0 0 0 0 0 0 0 0 0 0; +1 1 0 0 0 0 0 0 0 0 0; +1 1 1 0 0 0 0 0 0 0 0; +1 1 1 1 0 0 0 0 0 0 0; +1 1 1 1 1 0 0 0 0 0 0; +1 1 1 1 1 1 0 0 0 0 0; +1 1 1 1 1 1 1 0 0 0 0; +1 1 1 1 1 1 1 1 0 0 0; +1 1 1 1 1 1 1 1 1 0 0; +1 1 1 1 1 1 1 1 1 1 0]) + +groebner_walk(I, o_t, o_s, algorithm =:generic) #>1hr + diff --git a/experimental/GroebnerWalk/examples/knapsack/goodknap.jl b/experimental/GroebnerWalk/examples/knapsack/goodknap.jl new file mode 100644 index 000000000000..3ced896ba1fa --- /dev/null +++ b/experimental/GroebnerWalk/examples/knapsack/goodknap.jl @@ -0,0 +1,35 @@ +# "custom" integer knapsack problem + +using Oscar +R, (t, x1, x2, x3, x4, x5) = polynomial_ring(QQ, ["t","x1", "x2", "x3", "x4", "x5"]) + + + + + + +f1 = x1 - t^329221 +f2 = x2 - t^214884 +f3 = x3 - t^210568 +f4 = x4 - t^324307 +f5 = x5 - t^320945 + + +I = ideal([f1, f2, f3, f4,f5]) + +set_verbosity_level(:groebner_walk, 1) + +o_t = weight_ordering([1,0,0,0,0,0], degrevlex(R)) +o_s = weight_ordering([0,1,1,1,1,1], degrevlex(R)) +Ginit = groebner_basis(I, ordering = o_s, complete_reduction = true) +tg = @elapsed Gg = groebner_walk(I, o_t, o_s, algorithm =:generic) #60-70s + +ts = @elapsed Gs = groebner_walk(I, o_t, o_s, algorithm =:standard) + +tp = @elapsed Gp = groebner_walk(I, o_t, o_s, algorithm =:perturbed) +tb = @elapsed Gb = groebner_basis(I, ordering = o_t, complete_reduction = true) #140-160s + +tf = @elapsed Gf = groebner_basis(I, ordering = o_t, complete_reduction = true, algorithm =:fglm) #error: dimension of ideal must be zero! + +th = @elapsed Gh = groebner_basis(I, ordering = o_t, complete_reduction = true, algorithm =:hilbert) #throws another error + diff --git a/experimental/GroebnerWalk/examples/knapsack/knap.jl b/experimental/GroebnerWalk/examples/knapsack/knap.jl new file mode 100644 index 000000000000..6e1872d250d1 --- /dev/null +++ b/experimental/GroebnerWalk/examples/knapsack/knap.jl @@ -0,0 +1,20 @@ +using Oscar + +R, (t, x1, x2, x3, x4, x5, x6) = polynomial_ring(QQ, ["t","x1", "x2", "x3", "x4", "x5", "x6"]) + +f1 = x1 - t^53725 +f2 = x2 - t^61919 +f3 = x3 - t^64670 +f4 = x4 - t^69340 +f5 = x5 - t^78539 +f6 = x6 - t^95043 + +I = ideal([f1, f2, f3, f4, f5, f6]) + +set_verbosity_level(:groebner_walk, 1) + +o_t = weight_ordering([1,0,0,0,0,0,0], degrevlex(R)) +o_s = weight_ordering([0,1,1,1,1,1,1], degrevlex(R)) + +Gg = groebner_walk(I, o_t, o_s, algorithm =:generic) #300s + diff --git a/experimental/GroebnerWalk/examples/knapsack/smallknap.jl b/experimental/GroebnerWalk/examples/knapsack/smallknap.jl new file mode 100644 index 000000000000..9de58993e984 --- /dev/null +++ b/experimental/GroebnerWalk/examples/knapsack/smallknap.jl @@ -0,0 +1,18 @@ +#smallknap + +using Oscar +R, (t, x1, x2, x3) = polynomial_ring(QQ, ["t","x1", "x2", "x3"]) + +f1 = x1 - t^5 +f2 = x2 - t^12 +f3 = x3 - t^20 + +I = ideal([f1, f2, f3]) + +set_verbosity_level(:groebner_walk, 1) + +o_t = weight_ordering([1,0,0,0], degrevlex(R)) +o_s = weight_ordering([0,1,1,1], degrevlex(R)) + +Gg = groebner_walk(I, o_t, o_s, algorithm =:generic) + diff --git a/experimental/GroebnerWalk/examples/ku10.jl b/experimental/GroebnerWalk/examples/ku10.jl new file mode 100644 index 000000000000..33b51b108e36 --- /dev/null +++ b/experimental/GroebnerWalk/examples/ku10.jl @@ -0,0 +1,24 @@ +#ku +using Oscar + +R, (x1, x2, x3, x4, x5, x6, x7, x8, x9, x10) = polynomial_ring(QQ, ["x$index" for index in 1:10]) + +I = ideal([ + 5*x1*x2+ 5*x1+ 3*x2+ 55, + 7*x2*x3+ 9*x2+ 9*x3+ 19, + 3*x3*x4+ 6*x3+ 5*x4-4, + 6*x4*x5+ 6*x4+ 7*x5+ 118, + x5*x6+ 3*x5+ 9*x6+ 27, + 6*x6*x7+ 7*x6+x7+ 72, + 9*x7*x8+ 7*x7+x8+ 35, + 4*x8*x9+ 4*x8+ 6*x9+ 16, + 8*x9*x10+ 4*x9+ 3*x10-51, + 3*x1*x10-6*x1+x10+ 5 +]) + +set_verbosity_level(:groebner_walk, 1) + +t_s = @elapsed Gs = groebner_walk(I, lex(R), algorithm =:standard) +t_g = @elapsed Gg = groebner_walk(I, lex(R), algorithm =:generic) +t_p = @elapsed Gp = groebner_walk(I, lex(R), algorithm =:perturbed) + diff --git a/experimental/GroebnerWalk/examples/randomQQ.jl b/experimental/GroebnerWalk/examples/randomQQ.jl new file mode 100644 index 000000000000..ab907b3b0a89 --- /dev/null +++ b/experimental/GroebnerWalk/examples/randomQQ.jl @@ -0,0 +1,31 @@ +#random polynomial system in 4 variables over Q, computed in M2 +using Oscar +R, (a, b, c, d) = polynomial_ring(QQ, ["a", "b", "c", "d"]) + + +f1 = 15013d^3 + 5213*c*d^2 + 15893*c^2*d - 2837*b*d^2 - 10238*c^3 + 13538*b*c*d - 4578*a*d^2 + 5035*b*c^2 + +-9975*b^2*d - 9777*a*c*d - 2177*b^2*c + 9963*a*c^2 - 15569*a*b*d + 2513*b^3 - 3108*a*b*c - 9437*a^2*d - 5291*a*b^2 - +13120*a^2*c + 714*a^2*b + 13507*a^3 + +f2 = 2733*d^3 + 4541*c*d^2 + 4333*c^2*d - 6525*b*d^2 + 633*c^3 + 5406*b*c*d + 1762*a*d^2 + 11272*b*c^2 + 14879*b^2*d - +12968*a*c*d + 12580*b^2*c + 15050*a*c^2 + 15360*a*b*d - 2276*b^3 + 4773*a*b*c - 10499*a^2*d - 15025*a*b^2 - +9910*a^2*c - 3196*a^2*b + 1874*a^3 + +f3 = 11293*d^3 + 6929*c*d^2 - 13730*c^2*d - 7748*b*d^2 + 7572*c^3 - 9254*b*c*d + 15477*a*d^2 + 8782*b*c^2 + 7914*b^2*d + +5188*a*c*d + 15934*b^2*c - 10006*a*c^2 + 8830*a*b*d + 11842*b^3 + 7559*a*b*c + 11735*a^2*d - 7116*a*b^2 + +8771*a^2*c + 9251*a^2*b + 2099*a^3 + +I = ideal([f1, f2, f3]) + + +Gdegrevlex = groebner_basis(I, complete_reduction = true) +set_verbosity_level(:groebner_walk, 1) + +t_b = @elapsed Gb = groebner_basis(I, ordering = lex(R), complete_reduction = true) #8.6s + +t_s = @elapsed Gs = groebner_walk(I, lex(R), algorithm =:standard) #11.3, one conversion + +t_g = @elapsed Gg = groebner_walk(I, lex(R), algorithm =:generic) #14.4 + +t_p = @elapsed Gp = groebner_walk(I, lex(R), algorithm =:perturbed) #11.49 + diff --git a/experimental/GroebnerWalk/examples/toy_example.jl b/experimental/GroebnerWalk/examples/toy_example.jl new file mode 100644 index 000000000000..ccfd4b8c5224 --- /dev/null +++ b/experimental/GroebnerWalk/examples/toy_example.jl @@ -0,0 +1,10 @@ +using Oscar + +R,(x,y) = polynomial_ring(QQ, ["x","y"]) +I = ideal([y^4+ x^3-x^2+x,x^4]) + +set_verbosity_level(:groebner_walk, 1) + +groebner_walk(I) +groebner_walk(I; algorithm=:generic) + diff --git a/experimental/GroebnerWalk/examples/zero-dimensional/cyclic5.jl b/experimental/GroebnerWalk/examples/zero-dimensional/cyclic5.jl new file mode 100644 index 000000000000..f7ea458eb86d --- /dev/null +++ b/experimental/GroebnerWalk/examples/zero-dimensional/cyclic5.jl @@ -0,0 +1,20 @@ +using Oscar + +R, (a, b, c, d, x) = polynomial_ring(QQ, ["a", "b", "c", "d", "x"]) + +I = ideal([ + a + b + c + d + x, + a * b + b * c + c * d + d * x + x * a, + a * b * c + b * c * d + c * d * x + d * x * a + x * a * b, + a * b * c * d + b * c * d * x + c * d * x * a + d * x * a * b + x * a * b * c, + a * b * c * d * x - 1 +]) + +o_s = degrevlex(R) +o_t = lex(R) + +set_verbosity_level(:groebner_walk, 1) + +Gs = groebner_walk(I, lex(R), algorithm=:standard) +Gg = groebner_walk(I, lex(R), algorithm=:generic) + diff --git a/experimental/GroebnerWalk/examples/zero-dimensional/cyclic7.jl b/experimental/GroebnerWalk/examples/zero-dimensional/cyclic7.jl new file mode 100644 index 000000000000..8f134ec6017f --- /dev/null +++ b/experimental/GroebnerWalk/examples/zero-dimensional/cyclic7.jl @@ -0,0 +1,29 @@ +#cyclic7 + +using Oscar +R, (z0, z1, z2, z3, z4, z5, z6) = polynomial_ring(QQ, ["z$index" for index in 0:6 ]) + +I = ideal([ + z0 + z1 + z2 + z3 + z4 + z5 + z6, + + z0*z1 + z1*z2 + z2*z3 + z3*z4 + z4*z5 + z5*z6 + z6*z0, + + z0*z1*z2 + z1*z2*z3 + z2*z3*z4 + z3*z4*z5 + z4*z5*z6 + z5*z6*z0 + z6*z0*z1, + + z0*z1*z2*z3 + z1*z2*z3*z4 + z2*z3*z4*z5 + z3*z4*z5*z6 + z4*z5*z6*z0 + + z5*z6*z0*z1 + z6*z0*z1*z2, + + z0*z1*z2*z3*z4 + z1*z2*z3*z4*z5 + z2*z3*z4*z5*z6 + z3*z4*z5*z6*z0 + + z4*z5*z6*z0*z1 + z5*z6*z0*z1*z2 + z6*z0*z1*z2*z3, + + z0*z1*z2*z3*z4*z5 + z1*z2*z3*z4*z5*z6 + z2*z3*z4*z5*z6*z0 + z3*z4*z5*z6*z0*z1 + + z4*z5*z6*z0*z1*z2 + z5*z6*z0*z1*z2*z3 + z6*z0*z1*z2*z3*z4, + + z0*z1*z2*z3*z4*z5*z6 - 1 +]) + +set_verbosity_level(:groebner_walk, 1) + +Gs = groebner_walk(I, lex(R), algorithm =:standard) +Gg = groebner_walk(I, lex(R), algorithm =:generic) + diff --git a/experimental/GroebnerWalk/examples/zero-dimensional/finitefield.jl b/experimental/GroebnerWalk/examples/zero-dimensional/finitefield.jl new file mode 100644 index 000000000000..8397dd032a0a --- /dev/null +++ b/experimental/GroebnerWalk/examples/zero-dimensional/finitefield.jl @@ -0,0 +1,17 @@ +#= + Example over finite field, from `Groebner.jl` docstring line 1232 +=# +using Oscar +R, (x1, x2, x3, x4) = polynomial_ring(GF(32003), ["x1", "x2", "x3", "x4"]) + +J = ideal(R, [ + x1 + 2 * x2 + 2 * x3 + 2 * x4 - 1, + x1^2 + 2 * x2^2 + 2 * x3^2 + 2 * x4^2 - x1, + 2 * x1 * x2 + 2 * x2 * x3 + 2 * x3 * x4 - x2, + x2^2 + 2 * x1 * x3 + 2 * x2 * x4 - x3 +]) + +t_s = @elapsed Gs = groebner_walk(J, lex(R), algorithm=:standard) #4.11 +t_g = @elapsed Gg = groebner_walk(J, lex(R), algorithm=:generic) #0.8s +t_p = @elapsed Gp = groebner_walk(J, lex(R), algorithm=:perturbed) #0.33s + diff --git a/experimental/GroebnerWalk/examples/zero-dimensional/katsura5.jl b/experimental/GroebnerWalk/examples/zero-dimensional/katsura5.jl new file mode 100644 index 000000000000..f5e9b31847f3 --- /dev/null +++ b/experimental/GroebnerWalk/examples/zero-dimensional/katsura5.jl @@ -0,0 +1,23 @@ + +using Oscar + +R, (x,y,z,t,u,v) = QQ[:x,:y,:z,:t,:u,:v] + +I = ideal([ + 2*x^2+2*y^2+2*z^2+2*t^2+2*u^2+v^2-v, + x*y+y*z+2*z*t+2*t*u+2*u*v-u, + 2*x*z+2*y*t+2*z*u+u^2+2*t*v-t, + 2*x*t+2*y*u+2*t*u+2*z*v-z, + t^2+2*x*v+2*y*v+2*z*v-y, + 2*x+2*y+2*z+2*t+2*u+v-1 +]) + +o_s = degrevlex(R) +o_t = lex(R) + +set_verbosity_level(:groebner_walk, 0) + +t_s = @elapsed Gs = groebner_walk(I, lex(R), algorithm=:standard) +t_g = @elapsed Gg = groebner_walk(I, lex(R), algorithm=:generic) +t_b = @elapsed Gb = groebner_basis(I, ordering=lex(R)) + diff --git a/experimental/GroebnerWalk/results.csv b/experimental/GroebnerWalk/results.csv new file mode 100644 index 000000000000..2e9e2a2edc1a --- /dev/null +++ b/experimental/GroebnerWalk/results.csv @@ -0,0 +1,30 @@ +name,standard_walk,generic_walk,perturbed_walk,buchberger +simple,0.002634119,0.001080056,,4.670857142857143e-6 +simple,0.00262981,0.001087529,,4.621e-6 +simple,0.000623,0.000195084,,3.0837068965517243e-7 +cyclic5-QQ,0.01911975,0.110495292,,8.884516129032258e-7 +cyclic5-Fp,0.01191775,0.251762958,,8.766274509803922e-7 +simple,0.000620541,0.000199625,,3.2019369369369374e-7 +cyclic5-QQ,0.019888625,0.116872417,,8.674242424242424e-7 +cyclic5-Fp,0.01215325,0.268002959,,8.661276595744681e-7 +cyclic7-Fp,36.670678,,, +simple,0.000641,0.000202875,,3.0515983606557375e-7 +cyclic5-QQ,0.020051542,0.117901583,,8.731739130434783e-7 +cyclic5-Fp,0.013118917,0.261012125,,8.7975e-7 +cyclic6-QQ,0.356129125,5.922552333,,1.1916e-6 +cyclic6-Fp,0.130183917,11.620150625,,1.15e-6 +simple,0.000753417,0.000253458,,4.2939408866995077e-7 +simple,0.000753125,0.000237833,,4.0463111111111113e-7 +simple,1.183723958,0.269548625,,2.2375e-5 +agk4-QQ,10.771205958,55.705554,,527.126850292 +agk4-Fp,3.302066417,49.991833833,,0.095817208 +katsura6-QQ,1.003294291,,35.098429459, +simple,1.422464625,0.274453459,,2.475e-5 +agk4-QQ,11.076574625,55.562715875,,543.319831 +agk4-Fp,3.492181375,51.427781875,,0.094963834 +tran3.3-QQ,0.640738125,24.688099959,,2323.13229425 +tran3.3-Fp,0.254231917,14.402583584,,0.003486333 +newellp1-QQ,51.5992475,1529.494271334,,811.051479458 +newellp1-Fp,9.365612875,1140.743606416,,0.165012125 +simple,1.149968916,0.25867925,,2.5334e-5 +katsura6-QQ,0.813826709,35.404766583,, diff --git a/experimental/GroebnerWalk/src/GroebnerWalk.jl b/experimental/GroebnerWalk/src/GroebnerWalk.jl new file mode 100644 index 000000000000..5d68fb6a8ea8 --- /dev/null +++ b/experimental/GroebnerWalk/src/GroebnerWalk.jl @@ -0,0 +1,43 @@ +module GroebnerWalk +using Oscar + +import Oscar: + IdealGens, + MonomialOrdering, + ngens + weight_ordering, + ZZMatrix, + ZZRingElem + + +import Oscar.Orderings: + MatrixOrdering, + _support_indices + +include("markedGB.jl") +include("common.jl") + +include("special-ideals.jl") + +include("standard_walk.jl") +include("generic_walk.jl") +include("perturbed_walk.jl") + +export groebner_walk + +export newell_patch +export newell_patch_with_orderings + +function __init__() + add_verbosity_scope(:groebner_walk) +end + +end + +import .GroebnerWalk: + newell_patch, + newell_patch_with_orderings, + groebner_walk + +export groebner_walk + diff --git a/experimental/GroebnerWalk/src/common.jl b/experimental/GroebnerWalk/src/common.jl new file mode 100644 index 000000000000..089703a7d284 --- /dev/null +++ b/experimental/GroebnerWalk/src/common.jl @@ -0,0 +1,144 @@ + +@doc raw""" + groebner_walk( + I::MPolyIdeal, + target::MonomialOrdering = lex(base_ring(I)), + start::MonomialOrdering = default_ordering(base_ring(I)); + perturbation_degree = ngens(base_ring(I)), + algorithm::Symbol = :standard + ) + +Compute a reduced Groebner basis w.r.t. to the monomial ordering `target` by converting it +from a Groebner basis with respect to the ordering `start` using the Groebner Walk. + +# Arguments +- `I::MPolyIdeal`: ideal one wants to compute a Groebner basis for. +- `target::MonomialOrdering=:lex`: monomial ordering one wants to compute a Groebner basis for. +- `start::MonomialOrdering=:degrevlex`: monomial ordering to begin the conversion. +- `perturbationDegree::Int=2`: the perturbation degree for the perturbed Walk. +- `algorithm::Symbol=:standard`: strategy of the Groebner Walk. One can choose between: + - `standard`: Standard Walk [CLO05](@cite), + - `generic`: Generic Walk [FJLT07](@cite), + - `perturbed`: Perturbed Walk [AGK96](@cite). + +# Examples + +```jldoctest +julia> R,(x,y) = polynomial_ring(QQ, ["x","y"]); + +julia> I = ideal([y^4+ x^3-x^2+x,x^4]); + +julia> groebner_walk(I, lex(R)) +Gröbner basis with elements + 1 -> x + y^12 - y^8 + y^4 + 2 -> y^16 +with respect to the ordering + lex([x, y]) + +julia> groebner_walk(I, lex(R); algorithm=:generic) +Gröbner basis with elements + 1 -> y^16 + 2 -> x + y^12 - y^8 + y^4 +with respect to the ordering + lex([x, y]) + +julia> set_verbosity_level(:groebner_walk, 1); + +julia> groebner_walk(I, lex(R)) +Results for standard_walk +Crossed Cones in: +ZZRingElem[1, 1] +ZZRingElem[4, 3] +ZZRingElem[4, 1] +ZZRingElem[12, 1] +Cones crossed: 4 +Gröbner basis with elements + 1 -> x + y^12 - y^8 + y^4 + 2 -> y^16 +with respect to the ordering + lex([x, y]) + +``` +""" +function groebner_walk( + I::MPolyIdeal, + target::MonomialOrdering = lex(base_ring(I)), + start::MonomialOrdering = default_ordering(base_ring(I)); + perturbation_degree = ngens(base_ring(I)), # meaning, n=#gens(R) + algorithm::Symbol = :standard +) + if algorithm == :standard + walk = (x) -> standard_walk(x, target) + elseif algorithm == :generic + walk = (x) -> generic_walk(x, start, target) + elseif algorithm == :perturbed + walk = (x) -> perturbed_walk(x, start, target, perturbation_degree) + else + throw(NotImplementedError(:groebner_walk, algorithm)) + end + + Gb = groebner_basis(I; ordering=start, complete_reduction=true) + Gb = walk(Gb) + + return Oscar.IdealGens(Gb, target; isGB=true) +end + +@doc raw""" + is_same_groebner_cone(G::Oscar.IdealGens, T::MonomialOrdering) + +Checks whether the leading terms of G with respect to the matrix ordering given by T agree +with the leading terms of G with respect to the current ordering. +This means they are in the same cone of the Groebner fan. (cf. Lemma 2.2, Collart, Kalkbrener and Mall 1997) +""" +is_same_groebner_cone(G::Oscar.IdealGens, T::MonomialOrdering) = all(leading_term.(G; ordering=T) .== leading_term.(G; ordering=ordering(G))) + +# converts a vector wtemp by dividing the entries with gcd(wtemp). +@doc raw""" + convert_bounding_vector(w::Vector{T}) + +Scales a rational weight vector to have co-prime integer weights. +""" +convert_bounding_vector(w::Vector{T}) where {T<:Union{ZZRingElem, QQFieldElem}} = ZZ.(floor.(w//gcd(w))) + +# Creates a weight ordering for R with weight cw and representing matrix T +create_ordering(R::MPolyRing, cw::Vector{L}, T::Matrix{Int}) where {L<:Number} = weight_ordering(cw, matrix_ordering(R, T)) + +# interreduces the Groebner basis G. +# each element of G is replaced by its normal form w.r.t the other elements of G and the current monomial order +# TODO reference, docstring +@doc raw""" + autoreduce(G::Oscar.IdealGens) + +Replaces every element of G by the normal form with respect to the remaining elements of G and the current monomial ordering. +""" +function autoreduce(G::Oscar.IdealGens) + generators = collect(gens(G)) + + for i in 1:length(G) + generators[i] = reduce( + generators[i], generators[1:end .!= i]; ordering=G.ord, complete_reduction=true + ) + end + return Oscar.IdealGens(generators, G.ord; isGB=true) +end + +############################################# +# unspecific help functions +############################################# +#TODO docstring + +change_weight_vector(w::Vector{Int}, M::Matrix{Int}) = vcat(w', M[2:end, :]) +change_weight_vector(w::Vector{ZZRingElem}, M::ZZMatrix) = vcat(w', M[2:end, :]) + +insert_weight_vector(w::Vector{Int}, M::Matrix{Int}) = vcat(w', M[1:end-1, :]) +insert_weight_vector(w::Vector{ZZRingElem}, M::ZZMatrix) = vcat(w', M[1:end-1, :]) + +@doc raw""" + add_weight_vector(w::Vector{ZZRingElem}, M::ZZMatrix) + +Prepends the weight vector `w` as row to the matrix `M`. + +""" +add_weight_vector(w::Vector{Int}, M::Matrix{Int}) = vcat(w', M) +add_weight_vector(w::Vector{ZZRingElem}, M::ZZMatrix) = ZZMatrix(vcat(w), M) + diff --git a/experimental/GroebnerWalk/src/generic_walk.jl b/experimental/GroebnerWalk/src/generic_walk.jl new file mode 100644 index 000000000000..858922787831 --- /dev/null +++ b/experimental/GroebnerWalk/src/generic_walk.jl @@ -0,0 +1,246 @@ +@doc raw""" + generic_walk(G::Oscar.IdealGens, start::MonomialOrdering, target::MonomialOrdering) + +Compute a reduced Groebner basis w.r.t. to a monomial order by converting it using the Groebner Walk +using the algorithm proposed by [CKM97](@cite). + +# Arguments +- `G::Oscar.IdealGens`: Groebner basis of an ideal with respect to a starting monomial order. +- `target::MonomialOrdering`: monomial order one wants to compute a Groebner basis for. +- `start::MonomialOrdering`: monomial order to begin the conversion. +""" +function generic_walk(G::Oscar.IdealGens, start::MonomialOrdering, target::MonomialOrdering) + @vprintln :groebner_walk "Results for generic_walk" + @vprintln :groebner_walk "Facets crossed for: " + + Lm = leading_term.(G, ordering = start) + MG = MarkedGroebnerBasis(gens(G), Lm) + v = next_gamma(MG, ZZ.([0]), start, target) + + @v_do :groebner_walk steps = 0 + while !isempty(v) + @vprintln :groebner_walk v + MG = generic_step(MG, v, target) + v = next_gamma(MG, v, start, target) + + @v_do :groebner_walk steps += 1 + @vprintln :groebner_walk 2 G + @vprintln :groebner_walk 2 "=======" + end + + @vprint :groebner_walk "Cones crossed: " + @vprintln :groebner_walk steps + return gens(MG) +end + +@doc raw""" + generic_step(MG::MarkedGroebnerBasis, v::Vector{ZZRingElem}, ord::MonomialOrdering) + +Given the marked Gröbner basis `MG` and a facet normal vector `v`, compute the next marked Gröbner basis. +""" +function generic_step(MG::MarkedGroebnerBasis, v::Vector{ZZRingElem}, ord::MonomialOrdering) + facet_generators = facet_initials(MG, v) + H = groebner_basis( + ideal(facet_generators); ordering=ord, complete_reduction=true, algorithm=:buchberger + ) + + @vprint :groebner_walk 5 "Initials forms of facet: " + @vprintln :groebner_walk 5 facet_generators + + @vprint :groebner_walk 3 "GB of initial forms: " + @vprintln :groebner_walk 3 H + + H = MarkedGroebnerBasis(gens(H), leading_term.(H; ordering = ord)) + + lift_generic!(H, MG) + autoreduce!(H) + + return H +end + +@doc raw""" + difference_lead_tail(MG::MarkedGroebnerBasis) + +Computes $a - b$ for $a$ a leading exponent and $b$ in the tail of some $g\in MG$. +""" +function difference_lead_tail(MG::MarkedGroebnerBasis) + (G,Lm) = gens(MG), markings(MG) + lead_exp = leading_exponent_vector.(Lm) + + l_T = zip(lead_exp, exponents.(G)) + v = map(l_T) do (l,T) + Ref(l) .- T + end + + return [ZZ.(v)./ZZ(gcd(v)) for v in unique!(reduce(vcat, v)) if !iszero(v)] +end + +@doc raw""" + next_gamma( + MG::MarkedGroebnerBasis, w::Vector{ZZRingElem}, start::MonomialOrdering, target::MonomialOrdering + ) + +Given a marked Gröbner basis `MG`, a weight vector `w` and monomial orderings `start` and `target`, +returns the "next" facet normal, i.e. the bounding vector $v$ fulfilling $w maximum + max_deg = maximum(total_degree.(G)) + e = max_deg * m + 1 + @vprint :groebner_walk 5 "Upper degree bound: " + @vprintln :groebner_walk 5 e + + #w = sum(row * e^(p-i) for (i,row) in enumerate(rows)) + w = zeros(ZZRingElem, n) + d = 1 + for i in 1:p + w += row[end+1-i] * d + d *= e + end + + @vprint :groebner_walk 3 "Perturbed vector: " + @vprintln :groebner_walk 3 w + + return convert_bounding_vector(w) +end + diff --git a/experimental/GroebnerWalk/src/special-ideals.jl b/experimental/GroebnerWalk/src/special-ideals.jl new file mode 100644 index 000000000000..56f7b6226fe6 --- /dev/null +++ b/experimental/GroebnerWalk/src/special-ideals.jl @@ -0,0 +1,76 @@ + +@doc raw""" + newell_patch(k::Union{QQField, QQBarFieldElem}, n::Int=1) + +Let $n$ be an integer between 1 and 32. Returns the ideal corresponding to +the implicitization of the $n$-th bi-cubic patch describing +the Newell's teapot as a parametric surface. + +The specific generators for each patch have been taken from [Tra04](@cite). +""" +function newell_patch(k::Union{QQField, QQBarFieldElem}, n::Int=1) + return get_newell_patch_generators(n) |> ideal +end + +@doc raw""" + newell_patch(k::Field, n::Int=1) + +Let $n$ be an integer between 1 and 32. Returns the ideal corresponding to +the implicitization of the $n$-th bi-cubic patch describing +the Newell's teapot as a parametric surface. + +The specific generators for each patch have been taken from [Tra04](@cite). + +For fields $k\neq\mathbb{Q},\bar{\mathbb{Q}}$, this gives a variant of the ideal with +integer coefficients. +""" +function newell_patch(k::Field, n::Int=1) + F = get_newell_patch_generators(n) + + lcm_denom = [ + lcm(numerator.(coefficients(f))) for f in F + ] + integral_F = [ + change_coefficient_ring( + k, l*f + ) for (l, f) in zip(lcm_denom, F) + ] + + return ideal(integral_F) +end + +@doc raw""" + newell_patch_with_orderings(k::Field, n::Int=1) + +Let $n$ be an integer between 1 and 32. Returns the ideal corresponding to +the implicitization of the $n$-th bi-cubic patch describing +the Newell's teapot as a parametric surface. +Additionally returns suitable start and target orderings, e.g. for use with the Gröbner walk. + +The specific generators for each patch have been taken from [Tra04](@cite). + +For fields $k\neq\mathbb{Q},\bar{\mathbb{Q}}$, this gives a variant of the ideal with +integer coefficients. +""" +function newell_patch_with_orderings(k::Field, n::Int=1) + I = newell_patch(k, n) + R = base_ring(I) + o1 = matrix_ordering(R, [1 1 1 0 0; 0 0 0 1 1; 0 0 0 1 0; 1 1 0 0 0; 1 0 0 0 0]) + o2 = matrix_ordering(R, [0 0 0 1 1; 1 1 1 0 0; 1 1 0 0 0; 1 0 0 0 0; 0 0 0 1 0]) + return I, o2, o1 +end + +function get_newell_patch_generators(n::Int) + R, (x,y,z,u,v) = polynomial_ring(QQ, ["x","y","z","u","v"]) + + if n == 1 + return [ + -x + 7//5 - 231//125 * v^2 + 39//80 * u^2 - 1//5 * u^3 + 99//400 * u * v^2 - 1287//2000 * u^2 * v^2 + 33//125 * u^3 * v^2 - 3//16 * u + 56//125 * v^3 - 3//50 * u * v^3 + 39//250 * u^2 * v^3 - 8//125 * u^3 * v^3, + -y + 63//125 * v^2 - 294//125 * v + 56//125 * v^3 - 819//1000 * u^2 * v + 42//125 * u^3 * v - 3//50 * u * v^3 + 351//2000 * u^2 * v^2 + 39//250 * u^2 * v^3 - 9//125 * u^3 * v^2 - 8//125 * u^3 * v^3, + -z + 12//5 - 63//160 * u^2 + 63//160 * u + ] + else # TODO: Add the other patches + # TODO: Throw error + end +end + diff --git a/experimental/GroebnerWalk/src/standard_walk.jl b/experimental/GroebnerWalk/src/standard_walk.jl new file mode 100644 index 000000000000..b054931fbcb8 --- /dev/null +++ b/experimental/GroebnerWalk/src/standard_walk.jl @@ -0,0 +1,178 @@ +@doc raw""" + standard_walk(G::Oscar.IdealGens, start::MonomialOrdering, target::MonomialOrdering) + +Compute a reduced Groebner basis w.r.t. to a monomial order by converting it using the Groebner Walk +using the algorithm proposed by [CKM97](@cite). + +# Arguments +- `G::Oscar.IdealGens`: Groebner basis of an ideal with respect to a starting monomial order. +- `target::MonomialOrdering`: monomial order one wants to compute a Groebner basis for. +- `start::MonomialOrdering`: monomial order to begin the conversion. +""" +function standard_walk( + G::Oscar.IdealGens, + target::MonomialOrdering +) + start = ordering(G) + start_weight = canonical_matrix(start)[1, :] + target_weight = canonical_matrix(target)[1, :] + + return standard_walk(G, target, start_weight, target_weight) +end + +@doc raw""" + standard_walk(G::Oscar.IdealGens, target::MonomialOrdering, current_weight::Vector{ZZRingElem}, target_weight::Vector{ZZRingElem}) + +Compute a reduced Groebner basis w.r.t. to the monomial order `target` by converting it using the Groebner Walk +using the algorithm proposed by [CKM97](@cite). + +# Arguments +- `G::Oscar.IdealGens`: Groebner basis of an ideal with respect to a starting monomial order. +- `target::MonomialOrdering`: monomial order one wants to compute a Groebner basis for. +- `current_weight::Vector{ZZRingElem}`: weight vector representing the starting monomial order. +- `target_weight::Vector{ZZRingElem}`: weight vector representing the target weight. +""" +function standard_walk( + ::Type{Oscar.IdealGens}, + G::Oscar.IdealGens, + target::MonomialOrdering, + current_weight::Vector{ZZRingElem}, + target_weight::Vector{ZZRingElem}; +) + @vprintln :groebner_walk "Results for standard_walk" + @vprintln :groebner_walk "Crossed Cones in: " + + @v_do :groebner_walk steps = 0 + + while current_weight != target_weight + @vprintln :groebner_walk current_weight + G = standard_step(G, current_weight, target) + + current_weight = next_weight(G, current_weight, target_weight) + + @v_do :groebner_walk steps += 1 + @vprintln :groebner_walk 2 G + @vprintln :groebner_walk 2 "=======" + end + + @vprint :groebner_walk "Cones crossed: " + @vprintln :groebner_walk steps + + return G +end + +standard_walk( + G::Oscar.IdealGens, + target::MonomialOrdering, + current_weight::Vector{ZZRingElem}, + target_weight::Vector{ZZRingElem}; +) = gens(standard_walk(Oscar.IdealGens, G, target, current_weight, target_weight)) + +############################################################### +# The standard step is used for the strategies standard and perturbed. +############################################################### + +function standard_step(G::Oscar.IdealGens, w::Vector{ZZRingElem}, target::MonomialOrdering) + current_ordering = ordering(G) + next = weight_ordering(w, target) + + Gw = ideal(initial_forms(G, w)) + @vprint :groebner_walk 3 "GB of initial forms: " + @vprintln :groebner_walk 3 Gw + + H = groebner_basis(Gw; ordering=next, complete_reduction=true) + H = Oscar.groebner_lift(gens(H), next, gens(G), current_ordering) + + @vprint :groebner_walk 10 "Lifted GB of initial forms: " + @vprintln :groebner_walk 10 H + + return Oscar.IdealGens(H, next) +end +standard_step(G::Oscar.IdealGens, w::Vector{Int}, T::Matrix{Int}) = standard_step(G, ZZ.(w), create_ordering(base_ring(G), w, T)) + +@doc raw""" + initial_form(f::MPolyRingElem, w::Vector{ZZRingElem}) + +Returns the initial form of `f` with respect to the weight vector `w`. +""" +function initial_form(f::MPolyRingElem, w::Vector{ZZRingElem}) + R = parent(f) + + ctx = MPolyBuildCtx(R) + + E = exponents(f) + WE = dot.(Ref(w), Vector{ZZRingElem}.(E)) + maxw = -inf + + for (e, c, we) in zip(E, coefficients(f), WE) + if we > maxw + finish(ctx) + maxw = we + end + if we == maxw + push_term!(ctx, c, e) + end + end + + return finish(ctx) +end + +@doc raw""" + initial_forms(G::Oscar.IdealGens, w::Vector{ZZRingElem}) + +Returns the initial form of each element in `G` with respect to the weight vector `w`. +""" +initial_forms(G::Oscar.IdealGens, w::Vector{ZZRingElem}) = initial_form.(G, Ref(w)) +initial_forms(G::Oscar.IdealGens, w::Vector{Int}) = initial_form.(G, Ref(ZZ.(w))) + +@doc raw""" + next_weight(G::Oscar.IdealGens, current::Vector{ZZRingElem}, target::Vector{ZZRingElem}) + +Returns the point furthest along the line segment conv(current,target) still in the starting cone +as described in Algorithm 5.2 on pg. 437 of "Using algebraic geometry" (Cox, Little, O'Shea, 2005). + +# Arguments +- G, a reduced GB w.r.t the current monomial order +- current, a weight vector in the current Gröbner cone (corresponding to G) +- target a target vector in the Gröbner cone of the target monomial order +""" +function next_weight(G::Oscar.IdealGens, current::Vector{ZZRingElem}, target::Vector{ZZRingElem}) + V = bounding_vectors(G) + @vprint :groebner_walk 5 "Bounding vectors: " + @vprintln :groebner_walk 5 V + + C = dot.(Ref(current), V) + T = dot.(Ref(target), V) + + tmin = minimum(c//(c-t) for (c,t) in zip(C,T) if t<0; init=1) + + @vprint :groebner_walk 5 "Next rational weights: " + @vprintln :groebner_walk 5 (QQ.(current) + tmin * QQ.(target-current)) + + result = QQ.(current) + tmin * QQ.(target-current) + if !iszero(result) + return convert_bounding_vector(QQ.(current) + tmin * QQ.(target-current)) + else + return zeros(ZZRingElem, length(result)) + end +end + +@doc raw""" + bounding_vectors(I::Oscar.IdealGens) + +Returns a list of "bounding vectors" of a Gröbner basis of `I`, as pairs of +"exponent vector of leading monomial" and "exponent vector of tail monomial". +The bounding vectors form an H-description of the Gröbner cone. +(cf. p. 437 [CLO05](@cite)) +""" +function bounding_vectors(I::Oscar.IdealGens) + # TODO: Marked Gröbner basis + gens_by_terms = terms.(I; ordering=ordering(I)) + + v = map(Iterators.peel.(gens_by_terms)) do (lead,tail) + Ref(leading_exponent_vector(lead)) .- leading_exponent_vector.(tail) + end + + return unique!(reduce(vcat, v)) +end + diff --git a/experimental/GroebnerWalk/test/generic_walk.jl b/experimental/GroebnerWalk/test/generic_walk.jl new file mode 100644 index 000000000000..f793d04f03ea --- /dev/null +++ b/experimental/GroebnerWalk/test/generic_walk.jl @@ -0,0 +1,27 @@ +@testset "GenericWalk" begin + R1, (x,y,z) = polynomial_ring(QQ, ["x", "y", "z"]) + + I1 = ideal([x^2 + y*z, x*y + z^2]) + G1 = groebner_basis(I1) + + I2 = ideal([z^3 - z^2 - 4]) + G2 = groebner_basis(I2) + + I3 = ideal([x^6 - 1]) + G3 = groebner_basis(I3) + + GW1 = groebner_walk(I1; algorithm=:generic) + @test is_groebner_basis(GW1; ordering=lex(R1)) + #@test bounding_vectors(G1) == Vector{ZZRingElem}.([[1,1,-2],[2,-1,-1], [-1,2,-1]]) + #@test next_weight(G1, ZZ.([1,1,1]), ZZ.([1,0,0])) == ZZ.([1,1,1]) + + GW2 = groebner_walk(I2; algorithm=:generic) + @test is_groebner_basis(GW2; ordering=lex(R1)) + #@test bounding_vectors(G2) == ZZ.([[0,0,1],[0,0,3]]) + + GW3 = groebner_walk(I3; algorithm=:generic) + @test is_groebner_basis(GW3; ordering=lex(R1)) + #@test bounding_vectors(G3) == ZZ.([[6,0,0]]) + #@test next_weight(G3, ZZ.([1,1,1]), ZZ.([1,0,0])) == ZZ.([1,0,0]) +end + diff --git a/experimental/GroebnerWalk/test/standard_walk.jl b/experimental/GroebnerWalk/test/standard_walk.jl new file mode 100644 index 000000000000..590c2bf6556f --- /dev/null +++ b/experimental/GroebnerWalk/test/standard_walk.jl @@ -0,0 +1,24 @@ +@testset "StandardWalk" begin + R1, (x,y,z) = polynomial_ring(QQ, ["x", "y", "z"]) + + I1 = ideal([x^2 + y*z, x*y + z^2]) + G1 = groebner_basis(I1) + + I2 = ideal([z^3 - z^2 - 4]) + G2 = groebner_basis(I2) + + I3 = ideal([x^6 - 1]) + G3 = groebner_basis(I3) + + @test is_groebner_basis(groebner_walk(I1); ordering=lex(R1)) + #@test bounding_vectors(G1) == Vector{ZZRingElem}.([[1,1,-2],[2,-1,-1], [-1,2,-1]]) + #@test next_weight(G1, ZZ.([1,1,1]), ZZ.([1,0,0])) == ZZ.([1,1,1]) + + @test is_groebner_basis(groebner_walk(I2); ordering=lex(R1)) + #@test bounding_vectors(G2) == ZZ.([[0,0,1],[0,0,3]]) + + @test is_groebner_basis(groebner_walk(I3); ordering=lex(R1)) + #@test bounding_vectors(G3) == ZZ.([[6,0,0]]) + #@test next_weight(G3, ZZ.([1,1,1]), ZZ.([1,0,0])) == ZZ.([1,0,0]) +end + diff --git a/src/Rings/orderings.jl b/src/Rings/orderings.jl index 245dea1061cc..42fecc6e41ad 100644 --- a/src/Rings/orderings.jl +++ b/src/Rings/orderings.jl @@ -1232,7 +1232,7 @@ function matrix_ordering(v::AbstractVector{<:MPolyRingElem}, M::Union{Matrix{T}, end @doc raw""" - weight_ordering(W::Vector{Int}, ord::MonomialOrdering) -> MonomialOrdering + weight_ordering(W::Vector{<:IntegerUnion}, ord::MonomialOrdering) -> MonomialOrdering Given an integer vector `W` and a monomial ordering `ord` on a set of monomials in `length(W)` variables, return the monomial ordering `ord_W` on this set of monomials @@ -1280,7 +1280,7 @@ julia> canonical_matrix(o2) [0 0 1] ``` """ -function weight_ordering(w::Vector{Int}, o::MonomialOrdering) +function weight_ordering(w::Vector{<:IntegerUnion}, o::MonomialOrdering) i = _support_indices(o.o) m = ZZMatrix(1, length(w), w) return MonomialOrdering(base_ring(o), MatrixOrdering(i, m, false))*o